Abstract

Although the diversity of cortical interneuron electrical properties is well recognized, the number of distinct electrical types (e-types) is still a matter of debate. Recently, descriptions of interneuron variability were standardized by multiple laboratories on the basis of a subjective classification scheme as set out by the Petilla convention (Petilla Interneuron Nomenclature Group, PING). Here, we present a quantitative, statistical analysis of a database of nearly five hundred neurons manually annotated according to the PING nomenclature. For each cell, 38 features were extracted from responses to suprathreshold current stimuli and statistically analyzed to examine whether cortical interneurons subdivide into e-types. We showed that the partitioning into different e-types is indeed the major component of data variability. The analysis suggests refining the PING e-type classification to be hierarchical, whereby most variability is first captured within a coarse subpartition, and then subsequently divided into finer subpartitions. The coarse partition matches the well-known partitioning of interneurons into fast spiking and adapting cells. Finer subpartitions match the burst, continuous, and delayed subtypes. Additionally, our analysis enabled the ranking of features according to their ability to differentiate among e-types. We showed that our quantitative e-type assignment is more than 90% accurate and manages to catch several human errors.

Introduction

The diversity of cortical inhibitory interneurons is a well-accepted experimental phenomenon. This variability can be observed in various domains: morphological, that is, the anatomical structure of dendrites and axons (Ramon y Cajal 1904; Gibson et al. 1999; Cauli et al. 2000; McGarry et al. 2010), electrophysiological, that is, the firing patterns of the cells (McCormick et al. 1985; Connors and Gutnick 1990; Maccaferri and Lacaille 2003; McGarry et al. 2010), or the molecular, that is, the profile of gene expression. (Kawaguchi 1993; Cauli et al. 2000; Toledo-Rodriguez et al. 2004; Sugino et al. 2006). Different strategies have been employed in order to characterize this diversity, considering one or more of the aforementioned domain types (Gupta et al. 2000; Wang et al. 2002, 2004; Nowak et al. 2003; Toledo-Rodriguez et al. 2005; Halabisky et al. 2006; Ma et al. 2006; Helmstaedter et al. 2009a, 2009b; McGarry et al. 2010).

Recently, representatives from more than a dozen laboratories worldwide convened in Petilla de Aragon, Spain, for a discussion on the diversity of inhibitory interneurons (the “Petilla Interneuron Nomenclature Group” [PING]; Ascoli et al. 2008). The result of this convention was a unified nomenclature describing different forms of variability in interneurons.

When characterizing the diversity of interneuron electrical properties, many studies, including the PING, describe the diversity in terms of several distinct groups, or electrical type (e-types; sometimes also referred to as cell “classes”). In other words, a number of prototypical types of electrical properties (e.g., spiking response elicited by step currents) can be found from which each individual cell belonging to a class might vary to some extent, but considerably less than it varies from the prototype corresponding to a different e-type. Such a description is closely related to the notion of data clustering (Jain 2010). If the cells are indeed clustered into several distinct e-types rather than spread across a single continuum, the existence of distinct types can be demonstrated, and their number and nature can be objectively defined.

The caveat of the clustering approach is that one must define the “space” in which to check whether the data are clustered. This is especially problematic for electrophysiological data, as such a space is notoriously difficult to define. One main issue is the high-dimensional nature of the raw data. Spiking responses to particular stimuli typically require digitization at several kilohertz. If this were taken as the natural space, each sampling point would be a dimension in this space, and consequently there would be several thousand dimensions for even 1 s of voltage response recorded from a cell. This is clearly an unfeasible solution that is never performed in practice.

Alternatively, a set of descriptors, or features, of the voltage response to a particular stimulus can be extracted (e.g., the frequency of action potentials [APs], time of first AP following a particular current input, and AP width and height [Druckmann et al. 2007]). This set of features, together with clustering methods (Jain 2010), could then be used to explore whether the electrophysiological responses of the different cells fall into clusters corresponding to distinct types. The limitation of this approach is that the choice of features is somewhat arbitrary and that the transition from the raw voltage trace to the set of feature values involves information loss, as all the information not contained within the feature values is discarded. On the other hand, the manual PING classification may have its own limitations and biases. Therefore, it is only natural to study systematically, on the basis of a large set of features and quantitative statistical criteria, how robust the PING classification is and whether one can find additional structure within the data.

In this study, we extracted a large set of electrical features from the suprathreshold voltage response of several hundreds of cortical interneurons and analyzed this database of features using quantitative statistical approaches. The database of cortical interneurons was also manually annotated in accordance with the partitioning of types defined by the PING. In this way, we could both evaluate the explanatory power of the PING nomenclature for classifying e-types of cortical interneurons and examine which of the different features were most useful for capturing this partitioning of e-types.

We showed that the distribution of the electrophysiological feature values is clearly multimodal. These modes correspond quite well to the different types defined by the PING. Moreover, this coarse partition corresponds to the dimension in the data accounting for the highest variance. We showed that small subsets of features are sufficient to distinguish among these coarse types. Highly informative features are related both to the pattern of AP discharge (which is the basis of the naming convention in the PING) and to the shape of the individual APs. However, our analysis suggests a refinement for the PING classification, toward a more hierarchical view of the partitioning of cells into types. Notably, the aggregation of large sets of the data and its statistical analysis yielded a partitioning of the feature space into different “decision regions” by determining the e-type to which each point in space (representing a particular set of feature values) was most likely to belong. These decision regions can be used to quantitatively attribute a class to newly recorded cells, as well as to form the basis of comparisons with data from additional laboratories and different cell types.

Materials and Methods

Database of Electrophysiological Recordings

Data from in vitro intracellular recordings were collected as described in Toledo-Rodriguez et al. (2004, 2005) and aggregated into a large database of cortical neurons consisting of the response of 466 cortical interneurons from juvenile rats (age p 12–16). Features were calculated from responses to step currents of varying amplitude and length. In brief, different amplitude step currents were first employed to find the threshold current. Subsequent stimuli were then scaled according to this value. We analyze the responses to features extracted mainly from 2 amplitudes of step current injections: The first of amplitude equal to 150% threshold current and the second of amplitude equal to 300% threshold current (which we shall refer to as “standard” and “strong” current, respectively). Each step current had a duration of 2s. Approximately 50 traces were used for each neuron.

Electrical Features

We extracted a set of 38 features describing the cell's electrical response to step current pulse injection of varying amplitude. These features were based on a set of previously described electrical descriptors (Toledo-Rodriguez et al. 2004, 2005). Table 1 presents a brief description of each feature, a more detailed description of the features is found in the Supplementary Material. As features were defined in different units, prior to data analysis all features were transformed into z-scores, that is, the mean was subtracted and the feature value was normalized by the standard deviation (SD).

Table 1

Brief description of each feature

Feature # Feature description 
Drop in AP amplitude (amp.) from first to second spike (mV) 
AP amplitude change from first spike to steady-state (mV) 
AP 1 amplitude (mV) 
AP 1 width at half height (ms) 
AP 1 peak to trough time (ms) 
AP 1 peak to trough rate of change (mV/ms) 
AP 1 Fast AHP depth (mV) 
AP 2 amplitude (mV) 
AP 2 width at half height (ms) 
10 AP 2 peak to trough time (ms) 
11 AP 2 peak to trough rate of change (mV/ms) 
12 AP 2 Fast AHP depth (mV) 
13 Percent change in AP amplitude, first to second spike (%) 
14 Percent change in AP width at half height, first to second spike (%) 
15 Percent change in AP peak to trough rate of change, first to second spike (%) 
16 Percent change in AP fast AHP depth, first to second spike (%) 
17 Input resistance for steady-state current (Ohm) 
18 Average delay to AP 1 (ms) 
19 SD of delay to AP 1 (ms) 
20 Average delay to AP 2 (ms) 
21 SD of delay to AP 2 (ms) 
22 Average initial burst interval (ms) 
23 SD of average initial burst interval (ms) 
24 Average initial accommodation (%) 
25 Average steady-state accommodation (%) 
26 Rate of accommodation to steady-state (1/ms) 
27 Average accommodation at steady-state (%) 
28 Average rate of accommodation during steady-state 
29 Average inter-spike interval (ISI) coefficient of variation (CV) (unit less) 
30 Median of the distribution of ISIs (ms) 
31 Average change in ISIs during a burst (%) 
32 Average rate, strong stimulus (Hz) 
33 Average delay to AP 1, strong stimulus (ms) 
34 SD of delay to AP 1, strong stimulus (ms) 
35 Average delay to AP 2, strong stimulus (ms) 
36 SD of delay to AP 2, strong stimulus (ms) 
37 Average initial burst ISI, strong stimulus (ms) 
38 SD of average initial burst ISI, strong stimulus (ms) 
Feature # Feature description 
Drop in AP amplitude (amp.) from first to second spike (mV) 
AP amplitude change from first spike to steady-state (mV) 
AP 1 amplitude (mV) 
AP 1 width at half height (ms) 
AP 1 peak to trough time (ms) 
AP 1 peak to trough rate of change (mV/ms) 
AP 1 Fast AHP depth (mV) 
AP 2 amplitude (mV) 
AP 2 width at half height (ms) 
10 AP 2 peak to trough time (ms) 
11 AP 2 peak to trough rate of change (mV/ms) 
12 AP 2 Fast AHP depth (mV) 
13 Percent change in AP amplitude, first to second spike (%) 
14 Percent change in AP width at half height, first to second spike (%) 
15 Percent change in AP peak to trough rate of change, first to second spike (%) 
16 Percent change in AP fast AHP depth, first to second spike (%) 
17 Input resistance for steady-state current (Ohm) 
18 Average delay to AP 1 (ms) 
19 SD of delay to AP 1 (ms) 
20 Average delay to AP 2 (ms) 
21 SD of delay to AP 2 (ms) 
22 Average initial burst interval (ms) 
23 SD of average initial burst interval (ms) 
24 Average initial accommodation (%) 
25 Average steady-state accommodation (%) 
26 Rate of accommodation to steady-state (1/ms) 
27 Average accommodation at steady-state (%) 
28 Average rate of accommodation during steady-state 
29 Average inter-spike interval (ISI) coefficient of variation (CV) (unit less) 
30 Median of the distribution of ISIs (ms) 
31 Average change in ISIs during a burst (%) 
32 Average rate, strong stimulus (Hz) 
33 Average delay to AP 1, strong stimulus (ms) 
34 SD of delay to AP 1, strong stimulus (ms) 
35 Average delay to AP 2, strong stimulus (ms) 
36 SD of delay to AP 2, strong stimulus (ms) 
37 Average initial burst ISI, strong stimulus (ms) 
38 SD of average initial burst ISI, strong stimulus (ms) 

Multivariate Analysis

Following the extraction of the electrical descriptors, each cell was represented by a vector in m-dimensional space, where m = 38 is the number of features extracted (Table 1). As the data consisted of recordings from 466 cells, the data set was represented by an n × m matrix (n = 466, m = 38). Principal component analysis (PCA; Duda et al. 2001) was first performed to determine the prominent components of the variability in the data, by calculating the eigenvectors of the covariance matrix.

A linear Fisher discriminant analysis (Fisher 1936) was performed in order to determine the features that provide the greatest separation between the different types. Specifically, the discriminant functions (DFs) maximize the distance between the means of 2 groups, normalized by the scatter of each group. Formally, DFs are the eigenvectors of the following matrix: forumla where Sw, the within-group scatter matrix, was calculated by measuring the scatter of the feature values of cells belonging to each group following the subtraction from each feature of the group mean, and Sb is the scatter of the group means.

Feature subset selection was performed either by exhaustive search for up to 7 features or by a branch-and-bound algorithm (Narendra and Fukunaga 1977). In brief, branch-and-bound algorithms allow efficient evaluation of subsets of features if the criterion, in this case group separation, is monotonic in the features. Informally, monotonic functions are functions for which discarding additional features will only diminish the value of the criterion. These algorithms take advantage of the fact that if the criterion has already dropped below a certain value, throwing out more features is not going to improve it. Thus, when proceeding by eliminating (not adding) features, if a higher criterion value has been found for a different subset with a smaller number of features, then the lower-criterion subset need not be considered any more. Importantly, along with the lower-criterion subset, one can also reject without further consideration all the subsets derived from it by discarding additional features. Accordingly, many feature subsets do not need to be evaluated and the search becomes much more efficient than exhaustive search.

Decision Regions

In order to infer to which class a data point was most likely to belong, we needed to make some assumptions about the structure of the data. We assumed that each class could be treated as a multidimensional Gaussian, with different mean vectors and covariance matrices for each group. The parameters of each Gaussian (mean vector and covariance matrix) were estimated in standard fashion by maximum likelihood estimation (Duda et al. 2001): 

formula
 
formula

Given the Gaussian models and the estimated mean vector and covariance matrix, the probability for each point in space is defined by the following equation: 

formula
where |∑| indicates the determinant of the covariance matrix, sigma. This probability is calculated separately for each of the different Gaussians corresponding to the different classes. If all classes were equally common, every point in space would be assigned a class according to the Gaussian that it was most likely to arise from. When classes are not equally common, one must also factor in the relative contribution of each class to the population and each point is assigned a class according to the highest posterior probability: 
formula

The lines of equal posterior probability between groups are marked as the decision boundaries in Figure 4.

Testing for Normality

The classification analysis in the section above modeled each class as a 2-dimensional Gaussian. Accordingly, the closer the data match the Gaussian distribution, the more accurate the expected results. We tested for multivariate Gaussianity using Mardia's test for skewness and kurtosis (Mardia 1970).

Objective Nested Clustering Analysis

We performed a nested clustering analysis by repeatedly performing k-means clustering. Conceptually, this process is the methodic application of the process performed manually transitioning from Figures 2 and 3 (i.e., separating data sets and performing PCA on individual parts of the data set). Namely, beginning with the full data set, we calculated the first 10 PCs (that together capture 80% of the variance) and then performed k-means clustering to find 2 clusters within the data. These clusters are then split. For each of the clusters, we repeated the PCA with only the data from the individual cluster in order to reduce sensitivity to electrical class induced feature correlations. Then we ran separate k-means clustering on each of the clusters and split the cluster for which the splitting reduced the average Euclidean distance in PCA space between the group cluster center and each data point by the largest amount. For each of the k-means steps, we performed 500 repetitions of the k-means algorithm with different random seeds and selected the best cluster. We employed 10-fold cross validation to estimate the reliability of the results, that is, the process was repeated 10 times on 90% of the data each time. We note that we opted for this approach over more standard forms of divisive hierarchical clustering, such as (Macnaughton-Smith et al. 1964; Guénoche et al. 1991) due to the simplicity of k-means, its intuitive logic, and its widespread use in neuroscience and beyond. The average variance captured in the first 10 PCs at each of the splits is: (0.811, 0.825, 0.866, 0.869, 0.873, 0.916) respectively.

To check the consistency of the approach, we employed the Rand index (Rand 1971). In brief, the Rand index is an approach to measuring consistency between 2 sets of assignments of elements into classes. In our case, the 2 assignments are 2 repetitions of the nested clustering procedure, the classes are the e-types and the elements are the neurons. For every pair of neurons, one determines whether these 2 neurons are in the same class or not, for each of the 2 assignments. If the neurons are in the same class both in the first assignment and in the second, this pair is consistent across the assignments. Similarly, if the neurons are each in a different class both in the first assignment and in the second assignment, this pair is consistent as well. However, if the 2 points are in the same class in one assignment but not in the other, this is an inconsistency. Though very intuitive, the Rand index has one major disadvantage: It has no clear scale. Namely, though complete agreement between 2 assignments yields a score of 1, the expected value of the agreement between random assignments does not yield a specific score (such as zero). Therefore, we employed the adjusted Rand index that uses the generalized hypergeometric distribution as the control random model (Hubert and Arabie 1985).

Software

Data were queried from the database and imported into Matlab (Mathworks, Natick, MA, United States of America). Feature extraction was performed on all traces by custom software programmed in Matlab. Statistics of features and their normalization, as well as PCA were performed using built-in functions in Matlab. Multivariate analysis, including DF analysis and generation of decision regions, was programmed in Matlab.

Results

Figure 1 depicts the subjective classification of our database of cortical interneurons to e-types. Out of n = 466 cells recorded in our database, the majority (418) were subjectively classified in accordance with the Petilla (PING) nomenclature (Ascoli et al. 2008). For most of the e-types defined in this convention, we have a reasonable number of examples within the data set (Fig. 1). The PING convention delineates the responses into the steady-state part of the response, and the transient part of the response. In Figure 1, the steady-state responses comprise the rows (indicated by the following capital letters: fast spiking [FS], non-adapting [NA] non-fast spiking, adapting [AD], and irregular spiking, [IS]). The transient responses comprise the columns (indicated by lower case letters: burst [b], continuous [c], delayed [d], and stuttering [s]). In addition, we have the pyramidal cells (indicated by: Pyr). For 3 of the e-types (“burst-IS,” “continuous-IS,” and “delayed-non-adapting non-FS”) the database consists of less than 10 cells (gray traces). Accordingly, we did not treat these cells as separate types but rather labeled them collectively as the “other” cell type. Intrinsically bursting neurons were not considered because of the lack of reliable examples within the data considered for this study. A few cells (14) were marked as possibly belonging to 1 of 2 types and an unambiguous subjective classification was not found for a number of cells (31). Accordingly, these were treated as the “unspecified” type.

Figure 1.

Electrical cell types—subjective classification based on the PING convention. Color traces depict a single voltage trace in response to a depolarizing step current injection of intensity equal to 1.5 times the firing threshold intensity. One example is shown for all 10 interneuron electrical classes as well as one trace for pyramidal cells (black trace). The number of cells within the data set is shown for each group under the voltage trace. Gray traces are for cases containing less than 10 cells in the database and were thus pooled together into an “others” group.

Figure 1.

Electrical cell types—subjective classification based on the PING convention. Color traces depict a single voltage trace in response to a depolarizing step current injection of intensity equal to 1.5 times the firing threshold intensity. One example is shown for all 10 interneuron electrical classes as well as one trace for pyramidal cells (black trace). The number of cells within the data set is shown for each group under the voltage trace. Gray traces are for cases containing less than 10 cells in the database and were thus pooled together into an “others” group.

PCA of Feature Space

We began our unsupervised analysis by performing PCA based on the 38 features that we used for characterization of the electrical properties of each cell in our database (Table 1), using all the 466 neurons. Figure 2a shows the distribution of cells for different values of the first PC; the total distribution is depicted by the dotted black line, while the distribution for each of the subjectively classified e-types is marked by a corresponding colored line (AD, red, FS, blue, NA, green, Pyr, black, unspecified, orange, other, cyan, see legend). The total distribution function is clearly multimodal. The first 2 peaks correspond well to the FS cell type and to the AD e-type, respectively. The third peak is shared by Pyr neurons and adapting cells, and the remaining cell types are distributed fairly uniformly in this 1-dimensional projection upon the first PC.

Figure 2.

Cell classification based on PCA. (a) Distribution of cells as a function of the first PC score based on the 38 features used to characterize the cell's electrical properties (Table 1). Overall cell distribution is depicted by the dashed black line. Different subjectively classified e-types, according to PING convention, are marked by colored lines. Classes are marked with different colors as shown in inset (FS, fast spiking; AD, adapting; NA, non-adapting non-FS; Stut., stuttering; Pyr, pyramidal; Other, other; Unspec., unspecified). (b) Each cell is plotted as a point in the space of the first 2 PCs (PC1 and PC2). Point color is as in a, transient e-type is indicated by marker shape (c) Percentage of variance captured by each PC, PCs 1–10. (d) The weight of the contribution of each of the 38 features to the first (top) and second (bottom) PC. Features related to AP shape are marked in gray, whereas features related to AP firing pattern are marked in dark blue. Feature descriptions are provided in Table 1.

Figure 2.

Cell classification based on PCA. (a) Distribution of cells as a function of the first PC score based on the 38 features used to characterize the cell's electrical properties (Table 1). Overall cell distribution is depicted by the dashed black line. Different subjectively classified e-types, according to PING convention, are marked by colored lines. Classes are marked with different colors as shown in inset (FS, fast spiking; AD, adapting; NA, non-adapting non-FS; Stut., stuttering; Pyr, pyramidal; Other, other; Unspec., unspecified). (b) Each cell is plotted as a point in the space of the first 2 PCs (PC1 and PC2). Point color is as in a, transient e-type is indicated by marker shape (c) Percentage of variance captured by each PC, PCs 1–10. (d) The weight of the contribution of each of the 38 features to the first (top) and second (bottom) PC. Features related to AP shape are marked in gray, whereas features related to AP firing pattern are marked in dark blue. Feature descriptions are provided in Table 1.

In Figure 2b we plotted the cells in the space of the first 2 PCs. Each cell is represented as a point corresponding to its value in the first 2 PCs. Each point is colored according to the cell's subjective e-type classification, as in Figure 2a, and marker shape indicates transient e-type (circle—continuous, triangle—delayed, and square—bursting). As above, 2 visible clusters can be seen for the FS and adapting cells. Pyr cells occupy a diffuse region in the lower right part. Nearly all subjectively classified FS cells (blue) had negative values in the first PC, whereas adapting cells (red) had slightly negative or near-zero values and Pyr cells (black) obtained more positive values. Thus, the first PC is strongly related to the variability found in the subjective evaluation of cell types according to their steady-state electrical behavior (the rows of Fig. 1).

The fraction of the variability represented by each of the first 10 PCs is plotted in Figure 2c. Note that the first component is considerably larger than the subsequent ones. Thus, we confirm that the cell's steady-state e-type as defined by the PING nomenclature is indeed the major component of variability in the data. The data appear to be of a truly high-dimensional nature, with 10 PCs required to account for 80% of the variability in the data. In other words, there are numerous, separate factors that contribute to the electrical variability of cortical interneurons. Note that only the first 2 PCs are shown because of the limitations of plotting high-dimensional data.

The distribution of the weights of features associated with the first 2 PCs is depicted in Figure 2d. Features related to AP shape are marked in gray and features related to AP firing pattern in dark blue. Notably, the distribution is broad, with many of the features carrying substantial weight, both for features related to AP shape and to AP firing pattern [e.g., feature #4, width of the first AP at half height, and #30 median value of inter-spike-intervals (ISIs), corresponding to shape- and pattern-related features, respectively Table 1). Some features have only small weight [e.g., features #12 and #34 correspond to fast after-hyperpolariztion (AHP) of second AP, and SD of the delay to first AP in a strong stimulus] indicating that they do not contribute strongly to the main source of the variability in the data and consequently to the discrimination between the different steady-state neuron types.

We note that the different transient e-types do not stand out as separate clusters in Figure 2b. This could occur for 2 reasons: Either these different types are not in fact different types in the full high-dimensional feature space, or the particular projection chosen by the first 2 PCs causes these types to appear mixed. Why would the latter option be true? Judging by eye, the variability between the steady-state types is clearly large and immediately evident, whereas the differences between the transient behaviors are almost by definition more subtle. Moreover, the vast majority of cells come from the “continuous” transient behavior (see numbers in Fig. 1), further biasing the variance in the data toward differences between the steady-state e-types rather than those between transient types. Thus, PCA favors those directions in the high-dimensional space that correspond to separation between steady-state e-types. If these directions in feature space are different from those that separate between transient e-types, then the latter separation may be masked by the PCA (Supplementary Figure 1). We note that this is not a limitation of the PCA technique per se, but rather of using PCA (whose purpose is to find the directions in space that capture the most variance) to explore a different question, one of discrimination between e-types.

To test whether the differences in transient electrical behavior are masked by the differences in the steady-state electrical behaviors, we repeated PCA on data from a single steady-state e-type—the adapting interneurons, 158 neurons. The distribution of cells across the new first PC appears multimodal, though the 2 peaks overlap to some degree (Fig. 3a). The 2 peaks mostly correspond to the 2 transient electrical behaviors (burst transient and continuous) found for the adapting cells.

Figure 3.

PCA on adapting interneurons alone. (a) Spread of adapting interneuron feature values across the first PC. Black solid trace shows all adapting interneurons, dotted red trace shows cAD interneurons, and red solid line bAD interneurons. Note that the solid trace was artificially offset upwards to allow clearer viewing of all lines. (b) Weights assigned to different features composing the first PC. Features related to AP shape are shown in pink and features related to AP discharge pattern in dark red. Note that larger weights were assigned on average to the latter features (and compare also to Fig. 2d top). (c) Relative weight of the firing pattern features in different PCs. Weight is defined as the sum of the absolute value of feature weights, and is normalized by the different number of firing pattern versus AP shape features. Left and middle shows original PCs 1 and 2 (Fig. 2d) and right shows the PC plotted in Figure 3b.

Figure 3.

PCA on adapting interneurons alone. (a) Spread of adapting interneuron feature values across the first PC. Black solid trace shows all adapting interneurons, dotted red trace shows cAD interneurons, and red solid line bAD interneurons. Note that the solid trace was artificially offset upwards to allow clearer viewing of all lines. (b) Weights assigned to different features composing the first PC. Features related to AP shape are shown in pink and features related to AP discharge pattern in dark red. Note that larger weights were assigned on average to the latter features (and compare also to Fig. 2d top). (c) Relative weight of the firing pattern features in different PCs. Weight is defined as the sum of the absolute value of feature weights, and is normalized by the different number of firing pattern versus AP shape features. Left and middle shows original PCs 1 and 2 (Fig. 2d) and right shows the PC plotted in Figure 3b.

To understand why the difference between the transient types was not immediately apparent in the original PCA, we compared the weights assigned to the different features by PCA performed on all cells with the weight assigned by PCA performed separately on the adapting interneurons. The weights of the PCs were quite similar (Fig. 3b), yet the relative contribution of features related to AP discharge pattern, as opposed to AP shape, was larger in the separate PCA for adapting interneurons (Fig. 3c). This is important since the differences between burst-adapting (bAD) and continuous-adapting (cAD) are much more distinct in the AP discharge pattern. A one-way analysis of variance (ANOVA) shows highly significant group differences when based on the discharge features (n = 158, f = 23.61, P < 1e–5), yet the shape features show a much weaker, borderline non significant effect (n = 158, f = 3.68, P = 0.06). In summary, the high weights put on features related to AP shape in the original PC analysis (due to the distinct difference in shape between FS and AD interneurons and the large number of such neurons) masks the mostly discharge related differences between bAD and cAD interneurons.

DF Analysis

Having shown that the neurons' e-type is a major component of the variability in the data, we next explored the partitioning of neurons into different e-types by performing DF analysis (DFA, see Materials and Methods). In brief, DFA is a supervised data analysis method that finds the features that lead to the sharpest distinction between different classes. We stress that this method requires stating beforehand which neurons belong to which e-type and is in this sense biased. As noted previously, manual classification was performed for most cells by an expert based on the PING subjective criteria.

We performed DFA using the supervised manual classification. We use the 8 e-types that had a sufficient number of neurons (cFS, dFS, sFS, bNA, cNA, bAD, cAD, Pyr, 8 e-types, 387 neurons). In Figure 4a, we plotted the cells in the space of the first 2 DFs. The FS cells (blue, top right) are well separated from the adapting (red) cells, and even more distant from the Pyr (black) cells. Most adapting cells are well separated from the Pyr cells, but a portion are mixed among the Pyr cells (upper left portion of plot, red circles mixed with black circles) and appear to be more of a continuation of the Pyr cells than part of the main adapting e-type. The difference between the 3 e-types is highly significant, P < 1e–5 for all pairwise differences. The non-adapting non-FS (green) e-type is almost in full overlap with the adapting cells. The separation between the transient behavior (burst/continuous/delayed) cannot be distinguished at this resolution.

Figure 4.

DF analysis of electrical classes. (a) Each cell is plotted as a point in the space of the first 2 DFs. Marker color corresponds to subjectively characterized e-type, see legend. (b) Each feature in the discriminant analysis shown in a is plotted as a circle in the x,y location corresponding to the CV (mean over SD) and ratio of the within-group variance to the total variance, respectively (i.e., higher y-values indicate better separability of groups using this particular feature). Features related to AP shape are marked in gray and features related to AP firing pattern in black. For feature definition, see Table 1. (c) DF is performed on the fast-spiking and stuttering cells alone (the “blue” classes). Each fast-spiking cell is represented as a single point. Transient class is marked by different symbols (triangle—delayed, circle—continuous). (d) Each feature is plotted as in (b), but now with CV and variance ratio calculated for fast-spiking cells alone. Note the different y-axis scale in b and d.

Figure 4.

DF analysis of electrical classes. (a) Each cell is plotted as a point in the space of the first 2 DFs. Marker color corresponds to subjectively characterized e-type, see legend. (b) Each feature in the discriminant analysis shown in a is plotted as a circle in the x,y location corresponding to the CV (mean over SD) and ratio of the within-group variance to the total variance, respectively (i.e., higher y-values indicate better separability of groups using this particular feature). Features related to AP shape are marked in gray and features related to AP firing pattern in black. For feature definition, see Table 1. (c) DF is performed on the fast-spiking and stuttering cells alone (the “blue” classes). Each fast-spiking cell is represented as a single point. Transient class is marked by different symbols (triangle—delayed, circle—continuous). (d) Each feature is plotted as in (b), but now with CV and variance ratio calculated for fast-spiking cells alone. Note the different y-axis scale in b and d.

Figure 4b shows the contribution of different features to the separation between the e-types. In general, features that are well suited to separating any pair of 2 groups are those that take very different values if the neuron belongs to one group or another, but are relatively constant within each group. Formally, this is expressed as the ratio of the between-group variance to the within-group variance (see Materials and Methods). A second desirable property is that the variable not be very noisy, as expressed, for example, by its coefficient of variation (CV). Accordingly, each feature is plotted as a circle in the x,y position determined by the feature CV (x-axis) and the ratio of the between-group variance to the total variance (y-axis). Thus, features that offer good separability between the different types are found on the upper left part of the plot, corresponding to high between-to-within group variance and low CV. Indeed, a group of 8 features stands out on the top part of the plot. These are (in descending order of the ratio between group variance and total variance): (1) AP 1 width at half height; (2) average rate of strong stimulus; (3) AP 1 fast AHP depth; (4) average initial burst interval; (5) median of the distribution of ISIs; (6) AP 1 peak to trough rate of change; (7) SD of average initial burst interval; and (8) average ISI CV. Notably, the 2 main features that have been traditionally used to characterize FS cells—spike width and AP firing rate—are represented within this set.

We next attempted to distinguish between the transient electrical behaviors by DFA. As these can be obscured by the large differences between the steady-state electrical cell types as shown earlier, we performed discriminant analysis anew separately for each steady-state type. Figure 4c shows this analysis for the FS cells (the “blue” class—cFS, dFS, sFS; 116 neurons). Each FS cell is plotted as a point corresponding to the 2 first newly calculated DFs. Now the different transient types can be fairly well separated into their 3 subtypes. Figure 4d shows which features best discriminate between the different subtypes in this case. The important features in this case include: AP 2 peak to trough rate of change, AP 1 amplitude, average ISI CV, SD of delay from first to second spike, and AP 1 peak to trough time. These are mostly different from the features that separated between the steady-state firing types. Note, however, that even for the best feature, the between-group variance is only approximately 15% of the total variance of that feature (Fig. 4d, y-value of top circle), as opposed to the best feature for the steady-state electrical behavior discrimination shown in Figure 3b, where 50% of its variance was explained by the between-group variance (Fig. 4b, y-value of top circle). This indicates that the separation between the transient e-types is much weaker than the separation between the steady-state e-types. Similar analysis was performed for the difference between the different AD subtypes (Supplementary Figure 2).

Reliability of Discriminant Analysis

In order to test how well the e-types are separated, we employed multiple class separation in the space of the first 2 DFs. Namely, the mean and SD of a 2-dimensional Gaussian was estimated for each class separately. In addition, each point in space was assigned to the class whose Gaussian was most likely to generate it. Figure 5a shows the partition of the space into the points assigned to the 3 major types found above (black—Pyr, red—AD, blue—FS) the decision boundary is marked by respectively colored lines. The cells from these 3 types (335 neurons, encompassing all transient types) are shown atop this grid as full large dots. Testing for the hypothesis that each class arises from a 2-dimensional Gaussian (see Materials and Methods), we found that the null-hypothesis is not rejected (P > 0.05) for the FS and Pyr classes, whereas for the AD cells the hypothesis is rejected (P < 0.01) likely because of the existence of a number of extreme values.

Figure 5.

Classification by decision regions. The space of the first 2 DFs is divided into regions corresponding to the e-type that the data point would be most likely to have come from (see Materials and Methods): Adapting (red fine dots), FS (blue fine dots), or Pyr cells (black fine dots). Position of cells is shown by circles, colored according to e-type. Note that several circles lie far outside their decision area. (b) An example voltage response to a step current of 1.5 times firing threshold intensity is shown for the adapting (red, left) and FS (blue, right) electrical classes. (c). Three examples of cells manually classified as 1 e-type, yet located deep within the decision area of another e-type, are shown. The left column shows cells manually classified as adapting cells (red dots) yet located within the selection region of (and thus statistically classified as) FS cells. The right column shows cells manually classified as FS and statistically classified as adapting. Comparing the voltage responses to those found in b, the statistical analysis clearly gave the correct classification and the manual classification is erroneous.

Figure 5.

Classification by decision regions. The space of the first 2 DFs is divided into regions corresponding to the e-type that the data point would be most likely to have come from (see Materials and Methods): Adapting (red fine dots), FS (blue fine dots), or Pyr cells (black fine dots). Position of cells is shown by circles, colored according to e-type. Note that several circles lie far outside their decision area. (b) An example voltage response to a step current of 1.5 times firing threshold intensity is shown for the adapting (red, left) and FS (blue, right) electrical classes. (c). Three examples of cells manually classified as 1 e-type, yet located deep within the decision area of another e-type, are shown. The left column shows cells manually classified as adapting cells (red dots) yet located within the selection region of (and thus statistically classified as) FS cells. The right column shows cells manually classified as FS and statistically classified as adapting. Comparing the voltage responses to those found in b, the statistical analysis clearly gave the correct classification and the manual classification is erroneous.

Although some cells lie deep within the region belonging to a single e-type, they were subjectively classified as belonging to a different electric type. Four FS cells (blue circles) are found deep within the adapting cell (red) region and 5 adapting cells are found deep within the FS cell region. Figure 5c shows the firing response of 6 cells that were classified as belonging to a different e-type. As can be seen by the firing pattern, each of these cells actually belongs to the opposite e-type, as indeed they were assigned by the statistical classification of the decision boundaries. We note that such manual misclassification of 9 cells out of more than 300 corresponds only to a small percentage of the total number of cells and is an expected part of any human process applied to hundreds of examples, e.g., neuron tracing (Helmstaedter et al. 2011). Most importantly, the fact that mistakes in classifications can be caught alleviates at least in part the concern of using a supervised method like DFA, since it directly shows that the method does not simply recapitulate any arbitrary human labeling. Rather, by using the statistical power afforded by classification of a large number of cells, reliable measures of separation between the different e-types were determined despite conflicting examples (i.e., examples that reside in the portion of space occupied by an e-type different from that originally labeled).

Once such a decision region is found, it can be used to define the e-type of newly measured cells without the need for further subjective classification. In order to assess the accuracy of such an approach, we randomly excluded a portion of the data and used the remaining data to generate the decision region. We then checked the validity of the statistical classification by ten-fold cross validation (see Materials and Methods). Namely, we set aside data that were not used to construct the decision region, and estimated error rates by comparing the manual and statistical assignment on this validation set. The average rate of error was 6%, [mean error 5.9%± 2.3% (SD), 10-fold cross validation], meaning that we expect to correctly classify with no further human intervention about 94% of cells that will be measured in the future. We note that although FS interneurons are commonly distinguished in experiments from Pyr cells by the properties of their electrical response, adapting cells are considered to be difficult to distinguish from Pyr cells on the basis of their response to electrical stimuli alone. However, we find that this separation, although slightly more difficult, can still be achieved using our objective method with an 89% accuracy (mean error 89.1% ± 3.5%, 10-fold cross validation).

Between-Group and Within-Group Correlations

The importance of recognizing and describing e-types is especially relevant to examining the correlations between different biophysical properties. Correlations in the value of biophysical properties are often treated as suggestive of a tradeoff or balance between these different properties. However, one should be careful to distinguish between correlations that exist because of differences between group means and other, more straightforward correlations, which exist regardless of whether or not neurons belong to the same or a different class. For instance, we considered the relationship between the firing rate of a cell and the width of its AP for FS and AD cells (a total of 273 neurons). When pooling over the 2 cell types, we found a highly significant linear relationship for the data (Fig. 6, r2 = 0.422, P < 0.0001, confidence interval [CI] of slope [−72.1, −50.4]). Accordingly, one may be tempted to think that there is some biophysical mechanism that allows an increase in firing rate but causes spike width to be narrower. However, when we considered each of the 2 e-types separately, we found that this correlation was not significant, with the CI of the slope containing both negative and positive values. Thus, the correlation is almost exclusively due to the between-group effects, that is, the difference in the mean of the groups. Indeed, within each class, neurons with differences of almost half a millisecond in spike width (which is nearly the difference between the group means) exhibited similar firing rates. Thus, the nature and interpretation of the correlation between firing rate and AP width is more complex. In other words, an analysis blind to the e-type might have mistakenly assumed that the correlation exists within each group and not just due to differences in mean group values, a point that should guide its biophysical study.

Figure 6.

Electrical class and feature correlations. FS cells (blue) and adapting cells (red) are each represented by a colored dot located according to the value of their spike width (x-axis) and rate of AP discharge (y-axis). The overall correlation between AP width and firing rate for all cells is highly significant (black regression line, r2-value in black top right). However, when each group is considered separately, the regression between the 2 variables is found to be not significant (blue and red regression lines and r2-values in top right).

Figure 6.

Electrical class and feature correlations. FS cells (blue) and adapting cells (red) are each represented by a colored dot located according to the value of their spike width (x-axis) and rate of AP discharge (y-axis). The overall correlation between AP width and firing rate for all cells is highly significant (black regression line, r2-value in black top right). However, when each group is considered separately, the regression between the 2 variables is found to be not significant (blue and red regression lines and r2-values in top right).

Relative Importance of Features

Until now, we described analyses that had access to all the features. However, using a large number of features may complicate visualization and interpretation of data. Therefore, one may also consider using only smaller sets of features, thus raising the question: which are the most important features? To answer this question, we attempted to maximize the distinction between the types of electrical behavior as before, but here we used a more limited subset, each time of a different size. For small subsets one can search through all possible combinations of features, while for larger subsets, the large number of feature combinations requires more sophisticated searching, such as the use of branch-and-bound algorithms (Narendra and Fukunaga 1977, see Materials and Methods). Alternatively, there are numerous feature subset selection approaches that can be used also for large sets of features (for a review see Saeys et al. 2007). We performed an exhaustive search for feature subsets of up to size 4 and found that very accurate classifications could be performed on the basis of only a few features. This is perhaps to be expected since many of the features measure similar properties (Table 1).

We quantified the separability of the 3 steady-state e-types considered above (FS, AD, Pyr, 335 neurons) by the distance of the group means over the SD of the distributions (standardized mean difference, SMD). The larger SD of the 2 for each group was chosen as the normalization. In the case of the full feature sets, these differences were: FS-AD: 2.58 units of SD, AD-Pyr: 2.98 units. Considering subsets composed of only a single feature, we found that the most effective feature, AP 1 width at half-height, allows a rough classification (FS-AD 1.43, AD-Pyr 1.51, SMD). The best set of 3 features resulted in nicely separated groups (Fig. 7, FS-AD 2.21, AD-Pyr 2.16, SMD) and included the following features: AP 2 width at half-height, firing rate, AP 1 fast AHP depth. Notably, this set shares features that describe both the shape of the AP and the pattern of AP discharge. We obtained a classification nearly as accurate as that of the full feature set by using only 4 features: width of AP 2 at half-height, rate of AP discharge, AP 1 fast AHP depth, time to second AP from stimulus onset (FS-AD 2.33, AD-Pyr 2.51, SMD). Rerunning the analysis while allowing only features from 1 of the 2 sets (AP shape, discharge pattern) resulted in less accurate classification, which again highlights the usefulness of incorporating an analysis of the shape of the AP even though the electric types were originally defined (and are still named) according to their AP discharge pattern.

Figure 7.

Feature selection for steady-state e-type discrimination. A feature selection procedure was used to choose the features that best discriminate between FS interneurons (blue), adapting interneurons (red) and Pyr neurons (black). Plot shows the distribution of neuron feature values upon the first DF composed of the one best feature (top) and 3 best features (bottom). Note the sharper separation of groups in the bottom plot.

Figure 7.

Feature selection for steady-state e-type discrimination. A feature selection procedure was used to choose the features that best discriminate between FS interneurons (blue), adapting interneurons (red) and Pyr neurons (black). Plot shows the distribution of neuron feature values upon the first DF composed of the one best feature (top) and 3 best features (bottom). Note the sharper separation of groups in the bottom plot.

A Hierarchical Picture of Interneuron Electrical Variability

Taken together, these results strongly suggest a hierarchical picture for the classification of cortical interneuron electrical variability. Namely, the steady-state interneuron e-types (FS and AD cells) are easy to distinguish, account for a large portion of the variability while the transient electrical behavior types have smaller between-group differences (Figs 1–4 and the appropriate "Results" sections). To test this notion, we calculated the between-group distances for all e-type pairs (Fig. 8). For each of the e-types, we calculated the average of each feature value. This average defines a single point in feature space, that is, a space where each direction in is a given feature, yielding dimensionality equal to the number of features. Then, the distance between all pairs of points is computed. Since many of the features are correlated, we use the Mahalanobis distance (Mahalanobis 1936) rather than Euclidean distance. The distances found between the different e-types are indeed consistent with a hierarchical picture (Fig. 8a). First, the steady-state e-types have large distances between themselves (FS to AD). Second, the transient e-types within one steady-state response have smaller distances between themselves, that is, the cFS to dFS distance is smaller than that between FS and AD. Third, the distance between the same transient behavior but different steady-state behavior is larger than the distance between different transient behavior and same steady-state behavior, that is, the distance between cFS and cAD is larger than that between cFS and dFS. We note that the same properties also hold for distances calculated on the basis of groups directly taken from the original subjective classification.

Figure 8.

Hierarchical view of electrical variability of cortical interneurons. (a) Color plot of difference of class mean feature values between each pair of e-types. (b) Proposed schematic of the hierarchical structure of cortical interneuron electrical variability based on the hierarchical classification results. IS is marked by a question mark due to small representation in data.

Figure 8.

Hierarchical view of electrical variability of cortical interneurons. (a) Color plot of difference of class mean feature values between each pair of e-types. (b) Proposed schematic of the hierarchical structure of cortical interneuron electrical variability based on the hierarchical classification results. IS is marked by a question mark due to small representation in data.

To further test this question, we employed an unsupervised nested clustering approach that allows for changes in the correlations between different features (see Materials and Methods). In brief, considering the unambiguously labeled interneurons (404 neurons), the first split caused the majority of interneurons to break away from the Pyr neurons (Supplementary Figure 3). The following split caused the majority of FS cells to break away from the adapting cells. Next, the adapting cells broke into 2 groups, one of which contained most of the bAD interneurons. Next, the FS interneurons broke into a group that contained most of the delayed-FS type. The non-adapting non-FS cells tended to be distributed evenly between the different groups. Results were consistent within 10-fold cross validation (e.g., removing part of the data and rerunning the analysis, see Materials and Methods). We quantified consistency by the adjusted Rand index (Hubert and Arabie 1985, see Materials and Methods). Briefly, the adjusted Rand index compares 2 assignments of data points to partitions and yields an expected value of 1 if the 2 assignments are identical and 0 if the assignments are random. The adjusted Rand index found was highly significant (mean 0.68, SD = 0.03, P < 0.0001). Groups were then named by examining the subjective classification of each cell belonging to the group (which we stress was not used during the classification) and selecting the e-type that had the largest contribution. Though excluded from the full analysis because of the limited number of examples within our data set, we note that the IS neurons appear to branch from the adapting interneurons, with the aforementioned caveat of the small sample size.

In summary, the emerging hierarchical picture of interneuron e-types (Fig. 8b) appears well supported from the electrophysiological variability standpoint.

Discussion

In this study, we found that by quantifying the nomenclature described by the PING meeting (Ascoli et al. 2008) and applying several analysis approaches, a quantitative, statistical evaluation of the e-types characterized in this meeting could be readily obtained. Historically, the existence of multiple types of cortical interneurons (as well as hippocampal interneurons, e.g., Somogyi and Klausberger 2005) was suggested very early on using morphological features (Ramon y Cajal 1904). Since that time, cell classification has generated much interest, with different measures being used for neuron classification, including their firing properties, or e-types (McCormick et al. 1985; Connors and Gutnick 1990; Maccaferri and Lacaille 2003), their morphological structure (Ramon y Cajal 1904; Gibson et al. 1999; Cauli et al. 2000), and their genetic profile (Kawaguchi 1993; Cauli et al. 2000; Toledo-Rodriguez et al. 2004; Sugino et al. 2006). Regarding the classification into e-types, different approaches to what defines a class (Somogyi and Klausberger 2005; Ascoli et al. 2008; Helmstaedter et al. 2009a; McGarry et al. 2010) have led to a considerable amount of confusion, resulting in the often expressed feeling that such classification is inherently subjective, depending on whether one is a “lumper” or a “splitter.”

Importantly, the translation of the seemingly straightforward question “how many classes succinctly describe a certain dataset” into a quantitative framework is problematic. This is the case since it implies 2 contradictory goals: Describing the data well and having a minimal number of classes. The first, describing the data accurately, can be formalized as reducing the average distance between each individual data point and its associated group center (in our case, each data point corresponds to a neuron and group center to e-type). The optimal solution to this first goal is achieved by having one group for each data point, resulting in zero group center distance. The second goal, having a minimal amount of classes, is maximized by having one group for all the data points. Thus, these 2 aims are clearly in conflict and the tension between “lumping” and “splitting” is not merely a subjective one. We note that different approaches and heuristics can be employed in the attempt to find reasonable answers to the question of number of types (Duda et al. 2001) but the problem cannot be unequivocally resolved unless the exact generative model for the data is known, which is typically not the case.

We demonstrated in this study, which is based on a particularly large experimental database, that a partition into different major e-types is the foremost source of explainable variability in the data and that this partitioning is extremely salient statistically both in significance and in effect size. Notably, this partitioning can be found and confirmed with the simplest of statistical approaches (PCA, Gaussian decision theory). Thus, although one can argue that the natural tendency toward being a “lumper” or a “splitter” will affect one's view of the structure of interneuron electrical variability, the steady-state electrical properties should most certainly not be thought of as one continuum.

Supervised and Unsupervised Classification

Our simple statistical criteria were able to correctly predict the steady-state e-type of a cell (FS or AD) with an ∼95% accuracy. Like previous studies, we found that combining supervised and unsupervised classification techniques is a powerful approach (Guerra et al. 2011). We used the unsupervised techniques (PCA) to verify, in an exploratory unbiased way, that the separation into different e-types is indeed a salient feature of the data and not some preconceived notion that we have imposed on it. Following that, we used supervised techniques (DFA) to better understand the strongest features in terms of separating the e-types and to construct decision regions for classification. Although the criteria were learnt through analysis of the subjective classification, they were robust enough to uncover cases of manual misclassification (Fig. 5). Importantly, this means that we are now able to automatically classify a newly recorded cell without the need for further manual, subjective intervention. In addition, the entire process of transition from voltage traces to electrical features to classification is automated and requires very modest computational power; it can thus be readily extended to additional large databases of experimental recordings. Finally, we used the nested clustering approach (chosen on the basis of our finding that different features support different classifications) to confirm, in an unsupervised manner, the partitioning into different e-types. We note that the manual annotation used in this paper was based on a single expert. More robust labeling by taking the consensus opinion among multiple experts could have led perhaps to even stronger results.

Features

One of the main sources of controversy regarding the classification of interneuron e-types is the fact that there is no “natural space” in which the electrical activity of the cells should be described and their clustering (or lack thereof) examined, since the raw voltage traces' data are too high dimensional (including as many dimensions as sampling points). Thus, some form of dimensionality reduction is required. Typically, this takes the form of a set of features (e.g., firing rate, AP height, input resistance, etc.) that are used to represent the full data (Druckmann et al. 2007), entailing a loss of information.

The set of possible features can be quite large and any particular choice of features is partially arbitrary. Our choice in the current paper was to include a large set of features related to both the shape of individual APs (e.g., AP width) and the firing pattern (e.g., firing rate). Individual features were chosen with the aim of broadly characterizing these 2 types of information, and with the goal of including features used in previous studies. We note that large sets of features can be used, since as we demonstrated in this study, when a well-defined question is addressed (such as separation between labeled e-types), appropriate statistical tools for feature selection are available to assist in choosing minimal sets of these features in a principled fashion. Importantly, different features are informative for different questions, such as the separation between steady-state e-types versus transient e-types. One notable set of features that is likely to contain interesting information but was not included in this study in detail is the characterization of the firing rate of a neuron as a function of input current (F-I curves).

For the classification of the 3 major types of cortical neurons (FS interneurons, adapting interneurons, and Pyr amidal cells), we found that a highly accurate classification could be achieved even when only using a few features (in our hands: width of AP at half-height, rate of AP discharge, depth of AHP, time to second AP from stimulus onset). We note that the set of features found to be effective might vary for different data sets.

Use of Morphological, Electrophysiological, and Molecular Descriptors

In this study, we considered the electrophysiological properties of interneurons as a basis for their analysis and classification. Qualitatively different and equally important data regarding interneuron identity should also be considered, most notably morphological features of axons and dendrites, connectivity patterns, and molecular properties. The combination of different types of descriptors is crucial and can serve to both sharpen and validate the distinction between different neurons (Helmstaedter et al. 2009a, 2009b). We believe that the appropriate strategy is to first consider each type of descriptor separately and then combine them at the final stage of the analysis. Independently analyzing different types of descriptors allows their statistical power to be clearly examined. For instance, we show that adapting interneurons can be accurately distinguished from Pyr neurons based on their electrophysiological response alone, even if these responses appear quite similar to the naked eye. Notably, had we included a morphological attribute, such as maximal dendritic branch length, this major difference between interneurons and Pyr neurons might have masked the more subtle differences in their electrophysiological characteristics. Ultimately we believe the existence of neuronal types may be most convincingly shown if descriptors of different types (e.g., electrical, morphological, genetic) are clustered separately and then combined, thereby demonstrating that the types derived from the separate analyses can be related to each other as well as their associations (one-to-one, one-to-many, etc.).

A Hierarchy of Electrical Interneuron Types

In summary, we demonstrated that the PING nomenclature could be evaluated and backed up by quantitative statistical analysis. We found that the central aspects of the subjective classification described by the PING mostly matched the variability in our data and that it indeed can be parsimoniously explained as a partitioning of cortical cells into at least 3 major steady-state electrical activity types. However, the non-adapting non-FS class appears, at least in our data, to not be a distinct class (e.g., Fig. 1). Within these types further distinctions can be made according to transient properties. These distinctions are finer in terms of the variability in the data that they account for and are more difficult to quantify because of sampling biases in favor of the more common behaviors. We thus propose modifying the previous description of e-type partitioning through a 2-dimensional table (steady-state and transient behavior as proposed by the PING), in favor of a hierarchical view with the steady-state behavior being the more dominant factor (Fig. 8). Importantly, the ability to mark subsets of cells and specifically target them for data collection (Hempel et al. 2007) offers effective tools to offset sampling biases inherent to the under-representation of rare e-types and further examine these finer partitions.

Could the approach described above be extended for cell types beyond the neocortical interneurons considered here or to different brain areas? Since the analysis is largely automatic, adding different electrical properties to the feature set (e.g., those tied to different stimuli; Druckmann et al. 2011) to adapt the analysis to future data sets is rather straightforward. This is an important future direction that we will continue to address. We hope that the strength of large experimental databases shown in the framework provided above will provide incentive for the standardization of data collection and experimental procedures. Combining such data from different labs will more robustly test the approach presented in this study and mark new directions for extending the analysis.

Supplementary Material

Supplementary material can be found at: http://www.cercor.oxfordjournals.org/.

Funding

This work was supported by the EPFL fund for the Blue Brain Project, the Israel Science Foundation and by the Gatsby Charitable Foundation.

Notes

Conflict of Interest: None declared.

References

Ascoli
GA
Alonso-Nanclares
L
Anderson
SA
Barrionuevo
G
Benavides-Piccione
R
Burkhalter
A
Buzsaki
G
Cauli
B
Defelipe
J
Fairen
A
, et al.  . 
Petilla terminology: nomenclature of features of GABAergic interneurons of the cerebral cortex
Nat Rev Neurosci
 , 
2008
, vol. 
9
 (pg. 
557
-
568
)
Cauli
B
Porter
JT
Tsuzuki
K
Lambolez
B
Rossier
J
Quenet
B
Audinat
E
Classification of fusiform neocortical interneurons based on unsupervised clustering
Proc Natl Acad Sci USA
 , 
2000
, vol. 
97
 (pg. 
6144
-
6149
)
Connors
BW
Gutnick
MJ
Intrinsic firing patterns of diverse neocortical neurons
Trends Neurosci
 , 
1990
, vol. 
13
 (pg. 
99
-
104
)
Druckmann
S
Banitt
Y
Gidon
A
Schürmann
F
Markram
H
Segev
I
A novel multiple objective optimization framework for constraining conductance-based neuron models by experimental data
Front Neurosci
 , 
2007
, vol. 
1
 (pg. 
7
-
18
)
Druckmann
S
Berger
TK
Schürmann
F
Hill
S
Markram
H
Segev
I
Effective stimuli for constructing reliable neuron models
PLoS Comput Biol
 , 
2011
, vol. 
7
 pg. 
e1002133
 
Duda
RO
Hart
PE
Stork
DG
Pattern classification
 , 
2001
New York
Wiley
Fisher
RA
The use of multiple measurements in taxonomic problems
Ann Eugenic
 , 
1936
, vol. 
7
 (pg. 
179
-
188
)
Gibson
JR
Beierlein
M
Connors
BW
Two networks of electrically coupled inhibitory neurons in neocortex
Nature
 , 
1999
, vol. 
402
 (pg. 
75
-
79
)
Guénoche
A
Hansen
P
Jaumard
B
Efficient algorithms for divisive hierarchical clustering with the diameter criterion
J Classif
 , 
1991
, vol. 
8
 (pg. 
5
-
30
)
Guerra
L
McGarry
LM
Robles
V
Bielza
C
Larranaga
P
Yuste
R
Comparison between supervised and unsupervised classifications of neuronal cell types: A case study
Dev Neurobiol
 , 
2011
, vol. 
71
 (pg. 
71
-
82
)
Gupta
A
Wang
Y
Markram
H
Organizing principles for a diversity of GABAergic interneurons and synapses in the neocortex
Science
 , 
2000
, vol. 
287
 (pg. 
273
-
278
)
Halabisky
B
Shen
F
Huguenard
JR
Prince
DA
Electrophysiological classification of somatostatin-positive interneurons in mouse sensorimotor cortex
J Neurophysiol
 , 
2006
, vol. 
96
 (pg. 
834
-
845
)
Helmstaedter
M
Briggman
KL
Denk
W
High-accuracy neurite reconstruction for high-throughput neuroanatomy
Nat Neurosci
 , 
2011
, vol. 
14
 (pg. 
1081
-
1088
)
Helmstaedter
M
Sakmann
B
Feldmeyer
D
L2/3 interneuron groups defined by multiparameter analysis of axonal projection, dendritic geometry, and electrical excitability
Cereb Cortex
 , 
2009a
, vol. 
19
 (pg. 
951
-
962
)
Helmstaedter
M
Sakmann
B
Feldmeyer
D
The relation between dendritic geometry, electrical excitability, and axonal projections of L2/3 interneurons in rat barrel cortex
Cereb Cortex
 , 
2009b
, vol. 
19
 (pg. 
938
-
950
)
Hempel
CM
Sugino
K
Nelson
SB
A manual method for the purification of fluorescently labeled neurons from the mammalian brain
Nat Protoc
 , 
2007
, vol. 
2
 (pg. 
2924
-
2929
)
Hubert
L
Arabie
P
Comparing partitions
J Classif
 , 
1985
, vol. 
2
 (pg. 
193
-
218
)
Jain
AK
Data clustering: 50 years beyond K-means
Pattern Recogn Lett
 , 
2010
, vol. 
31
 (pg. 
651
-
666
)
Kawaguchi
Y
Physiological, morphological, and histochemical characterization of three classes of interneurons in rat neostriatum
J Neurosci
 , 
1993
, vol. 
13
 (pg. 
4908
-
4923
)
Ma
Y
Hu
H
Berrebi
AS
Mathers
PH
Agmon
A
Distinct subtypes of somatostatin-containing neocortical interneurons revealed in transgenic mice
J Neurosci
 , 
2006
, vol. 
26
 (pg. 
5069
-
5082
)
Maccaferri
G
Lacaille
JC
Interneuron diversity series: Hippocampal interneuron classifications—making things as simple as possible, not simpler
Trends Neurosci
 , 
2003
, vol. 
26
 (pg. 
564
-
571
)
Macnaughton-Smith
P
Williams
WT
Dale
MB
Mockett
LG
Dissimilarity analysis: A new technic of hierarchical subdivision
Nature
 , 
1964
, vol. 
202
 (pg. 
1034
-
1035
)
Mahalanobis
PC
On the generalised distance in statistics
Proc Natl Inst Sci
 , 
1936
, vol. 
2
 (pg. 
49
-
55
)
Mardia
KV
Measures of multivariate skewness and kurtosis with applications
Biometrika
 , 
1970
, vol. 
57
 pg. 
519
 
McCormick
DA
Connors
BW
Lighthall
JW
Prince
DA
Comparative electrophysiology of pyramidal and sparsely spiny stellate neurons of the neocortex
J Neurophysiol
 , 
1985
, vol. 
54
 (pg. 
782
-
806
)
McGarry
LM
Packer
AM
Fino
E
Nikolenko
V
Sippy
T
Yuste
R
Quantitative classification of somatostatin-positive neocortical interneurons identifies three interneuron subtypes
Front Neural Circuits
 , 
2010
, vol. 
4
 pg. 
12
 
Narendra
PM
Fukunaga
K
A branch and bound algorithm for feature subset selection
IEEE Trans Comput
 , 
1977
, vol. 
100
 (pg. 
917
-
922
)
Nowak
LG
Azouz
R
Sanchez-Vives
MV
Gray
CM
McCormick
DA
Electrophysiological classes of cat primary visual cortical neurons in vivo as revealed by quantitative analyses
J Neurophysiol
 , 
2003
, vol. 
89
 (pg. 
1541
-
1566
)
Ramon y Cajal
S
Textura del sistema nervioso del hombre y de los vertebrados
 , 
1904
Madrid
Imprenta N. Moya
Rand
WM
Objective criteria for evaluation of clustering methods
J Am Stat Assoc
 , 
1971
, vol. 
66
 (pg. 
846
-
853
)
Saeys
Y
Inza
I
Larranaga
P
A review of feature selection techniques in bioinformatics
Bioinformatics
 , 
2007
, vol. 
23
 (pg. 
2507
-
2517
)
Somogyi
P
Klausberger
T
Defined types of cortical interneurone structure space and spike timing in the hippocampus
J Physiol
 , 
2005
, vol. 
562
 (pg. 
9
-
26
)
Sugino
K
Hempel
CM
Miller
MN
Hattox
AM
Shapiro
P
Wu
C
Huang
ZJ
Nelson
SB
Molecular taxonomy of major neuronal classes in the adult mouse forebrain
Nat Neurosci
 , 
2006
, vol. 
9
 (pg. 
99
-
107
)
Toledo-Rodriguez
M
Blumenfeld
B
Wu
C
Luo
J
Attali
B
Goodman
P
Markram
H
Correlation maps allow neuronal electrical properties to be predicted from single-cell gene expression profiles in rat neocortex
Cerebral Cortex
 , 
2004
, vol. 
14
 (pg. 
1310
-
1327
)
Toledo-Rodriguez
M
Goodman
P
Illic
M
Wu
C
Markram
H
Neuropeptide and calcium-binding protein gene expression profiles predict neuronal anatomical type in the juvenile rat
J Physiol
 , 
2005
, vol. 
567
 (pg. 
401
-
413
)
Wang
Y
Gupta
A
Toledo-Rodriguez
M
Wu
CZ
Markram
H
Anatomical, physiological, molecular and circuit properties of nest basket cells in the developing somatosensory cortex
Cereb Cortex
 , 
2002
, vol. 
12
 (pg. 
395
-
410
)
Wang
Y
Toledo-Rodriguez
M
Gupta
A
Wu
C
Silberberg
G
Luo
J
Markram
H
Anatomical, physiological and molecular properties of Martinotti cells in the somatosensory cortex of the juvenile rat
J Physiol
 , 
2004
, vol. 
561
 (pg. 
65
-
90
)