|Home | About | Journals | Submit | Contact Us | Français|
Deciphering the circuitry of the neocortex requires knowledge of its components, making a systematic classification of neocortical neurons necessary. GABAergic interneurons contribute most of the morphological, electrophysiological and molecular diversity of the cortex, yet interneuron subtypes are still not well defined. To quantitatively identify classes of interneurons, 59 GFP-positive interneurons from a somatostatin-positive mouse line were characterized by whole-cell recordings and anatomical reconstructions. For each neuron, we measured a series of physiological and morphological variables and analyzed these data using unsupervised classification methods. PCA and cluster analysis of morphological variables revealed three groups of cells: one comprised of Martinotti cells, and two other groups of interneurons with short asymmetric axons targeting layers 2/3 and bending medially. PCA and cluster analysis of electrophysiological variables also revealed the existence of these three groups of neurons, particularly with respect to action potential time course. These different morphological and electrophysiological characteristics could make each of these three interneuron subtypes particularly suited for a different function within the cortical circuit.
Despite comprising a minority of all neocortical neurons, GABAergic interneurons appear to play an important circuit role by responding to dynamic changes in excitation, perhaps to keep the circuit responsive over a wide range of inputs, synchronize activity, control timing of pyramidal cell firing, or suppress runaway excitation (McBain and Fisahn, 2001; Pouille and Scanziani, 2001; Markram et al., 2004; Trevelyan et al., 2006; Klausberger and Somogyi, 2008). In support of these roles, GABAergic neurons have been implicated in several pathologies including epilepsy (Cossart et al., 2001; Cobos et al., 2005; Freund and Katona, 2007; Trevelyan et al., 2006, 2007), autism (Tabuchi et al., 2007), Rett syndrome (Dani et al., 2005), anxiety disorder (Powell et al., 2003; Freund and Katona, 2007), Tourette syndrome (Kalanithi et al., 2005) and schizophrenia (Lewis et al., 2005; Gonzalez-Burgos and Lewis, 2008).
Interneurons exhibit very diverse morphology, electrophysiology, molecular content and post synaptic targets (Cauli et al., 1997; Markram et al., 2004; Yuste, 2005). Classification of neocortical interneurons is a crucial step in understanding cortical circuits as each subtype of interneuron likely has a different function. In the past, classifications have been based on one or a combination of descriptors, often using qualitative criteria which were not standardized. Therefore, interneuron classifications often differ and it is unclear which set of descriptors are the most relevant to determine a neuronal class and, more generally, how many classes of interneurons actually exist. As a first step, several groups (Cauli et al., 2000; Karube et al., 2004; Dumitriu et al., 2006; Ma et al., 2006; Helmstaedter et al., 2008; Karagiannis et al., 2009) have used quantitative methods, with unsupervised clustering algorithms to seek an objective classification of interneurons. In addition, to facilitate the classification of interneurons, a standardized nomenclature for morphological, physiological and molecular features was agreed upon by the Petilla Interneuron Nomenclature Group (Ascoli et al., 2008).
Among GABAergic interneurons, those expressing the neuropeptide somatostatin (SOM) are particularly heterogeneous in their molecular, morphological and electrophysiological features (Sabo and Sceniak, 2006). To better understand this variability, we took advantage of a transgenic mouse line where SOM-positive neurons are labeled with GFP (Oliva et al., 2000) to explore if subtypes of SOM interneurons could be distinguished objectively. 59 SOM-positive interneurons from mouse primary somatosensory, frontal and visual cortex were quantitatively characterized using their morphological properties measured from reconstructions of biocytin-filled cells and their intrinsic electrophysiological properties, measured from whole cell recordings. Unsupervised cluster analysis revealed a group comprised of the well-known Martinotti cells, as well as two other sub-groups of SOM-positive neurons that, to our knowledge, have not yet been described. The electrophysiological classification agreed well with, and thus confirmed, the morphological classification. Furthermore, significant differences in axon morphology and firing properties between the three groups suggest that the two novel subtypes and Martinotti cells could serve different functional roles in the cortical circuit.
Acute brain slices were prepared from GIN mice (Oliva et al., 2000), with an average age postnatal 14days (range P10–18). Mice were quickly decapitated, the skin and skull were cut and the brain was removed and then immediately placed in cold sucrose cutting solution (222mM sucrose, 2.6mM KCl, 27mM NaHCO3, 1.5mM NaH2PO4, 2mM CaCl2, 2mM MgSO4, bubbled with 95% 02, 5%CO2). Slices 300-μm thick were cut using a Vibratome and transferred to a chamber at room temperature with oxygenated ACSF (126mM NaCl, 3mM KCl, 2mM MgSO4, 2mM CaCl2, 1mM NaH2PO4, 26mM NaHCO3 and 10mM glucose, bubbled with 95% 02, 5%CO2).
GIN mice reliably label a subset of SOM-positive cells with expression remaining consistent for more than five generations (Oliva et al., 2000). Two recent studies have repeated immunohistochemical analysis of this original paper, demonstrating that GFP expressing cells are indeed SOM-positive. One group reported co-localization in 95.9% of layer 2/3 cells, 100% of layer 4 cells and 97.4% of layer 5/6 cells (Ma et al., 2006) and the other reported co-localization of 99% for cells in all layers (Halabisky et al., 2006). In this study, the percentage of SOM+cells that were GFP labeled was 34.8% in layer 2/3, 26.5% in layer 4 and 10.8% in layer 5/6. To our knowledge no one line of transgenic mice labels all SOM+cells and the most often used GIN line labels a higher percentage of SOM+cells than other lines (Ma et al., 2006).
Slices were placed in a recording chamber at room temperature with flowing oxygenated ACSF. Pipettes of 3–7MΩ resistance were pulled from borosilicate glass. SOM positive cells were identified by expression of GFP. Whole cell recordings were taken in current clamp mode. Only cells with healthy resting membrane potential (between −55 and −80mV) were selected for recording.
Ninteen variables were measured for each neuron by analysis of the recordings in MATLAB. See Table Table11 for descriptions.
Neurons were filled with biocytin by a patch pipette. The slices were kept overnight in 4% paraformaldehyde in 0.1M phosphate buffer (PB) at 4°C. The slices were then rinsed three times for 5min per rinse on a shaker in 0.1M PB. They were placed in 30% sucrose mixture (30g sucrose dissolved in 50ml ddH20 and 50ml 0.24M PB per 100ml) for 2h and then frozen on dry ice in tissue freezing medium. The slices were kept overnight in a −80°C freezer. The slices were defrosted and the tissue freezing medium was removed by three 20-min rinses in 0.1M PB while on a shaker. The slices were kept in 1% hydrogen peroxide in 0.1M PB for 30min on the shaker to pretreat the tissue, then were rinsed twice in 0.02M potassium phosphate saline (KPBS) for 20min on the shaker. The slices were then kept overnight on the shaker in Avidin-Biotin-Peroxidase Complex. The slices were then rinsed three times in 0.02M KPBS for 20min each on the shaker. Each slice was then placed in DAB (0.7mg/ml 3,3′-diaminobenzidine, 0.2mg/ml urea hydrogen peroxide, 0.0 6M Tris buffer in 0.02M KPBS) until the slice turned light brown then immediately transferred to 0.02M KPBS and transferred again to fresh 0.02M KPBS after a few minutes. The stained slices were rinsed a final time in 0.02M KPBS for 20min on a shaker. Each slice was observed under a light microscope and then mounted onto a slide using crystal mount.
Successfully filled and stained neurons were then reconstructed using Neurolucida software (MicroBrightField). The neurons were viewed with 100× oil objective on an Olympus IX71 inverted light microscope or an Olympus BX51 upright light microscope. The Neurolucida program projects the microscope image onto a computer drawing tablet. The neuron's processes were traced manually while the program recorded the coordinates of the tracing to create a digital three dimensional reconstruction. The x and y axis form the horizontal plane of the slice, while the z axis is the depth. The user defined an initial reference point for each tracing. The z coordinate was then determined by adjustment of the focus. In addition to the neuron, the pia and white matter were drawn.
The Neurolucida Explorer program was used to measure 67 morphological variables of the reconstruction. See Table Table22 for descriptions.
Shrinkage percentage in the z-axis was on average 48.97±1.95%. No shrinkage correction was applied.
The dataset of neurons could now be represented as a matrix, with each row corresponding to a neuron and each column corresponding to a variable. The data was then standardized, where the standardized values of each variable are computed as the difference between the actual value and the mean divided by the standard deviation. PCA reduces the dimensionality of the dataset, while preserving as much of the variance as possible by mapping the coordinate system defined by the old variables to a new lower dimensional coordinate system defined by principal components that are orthogonal, thus uncorrelated, linear combinations of the original variables. As many of the original variables are related, the principal components produced by PCA reduce this redundancy. By PCA, a new space is generated onto which the dataset can be projected and classified by cluster analysis.
The principal components were calculated using the factor analysis function in STATISTICA (StatSoft). The principal components are found by an eigenvalue decomposition of the correlation matrix of the standardized data. The eigenvectors of the correlation matrix are the principal components (PCs) and the corresponding eigenvalues give the variance preserved by each principal component. The number of PCs maintaining a significant amount of variance will always be less than the number of original variables, which is why PCA reduces the dimensionality of the dataset (Jambu, 1991).
Several criteria were used to decide how many principal components to retain for cluster analysis. First, only principal components with eigenvalues greater than one are considered. The original variables have an eigenvalue of one, as the data is standardized, thus any principal component with an eigenvalue greater than one describes more of the data's variance than an original variable. Second, the “scree test” was used (Cattell, 1966). A plot of the eigenvalues of the principal components is examined for when the decrease in eigenvalue plateaus. In all cases either 2 or 3 principal components were retained.
Cluster Analysis is then performed on each dataset in the principal component vector space using Clustan (Clustan Ltd.). Hierarchical clustering using Euclidean distance squared as the multi-dimensional distance metric and Ward's method as the linkage rule was performed. Ward's method minimizes the distance squared between any two clusters that can be formed at each step of clustering (Ward, 1963). Alternatively, cluster analysis was performed on the standardized data before PCA. The results obtained from this method did not dramatically differ from the results obtained with PCA. The only differences were reordering of cells within clusters and movement of cells between neighboring clusters contained in Group 1.
Three statistical tests were done to find the clustering levels which have significance at the 5% level. The Best Cut Test does an analysis of variance on the fusion values (distance at which two clusters join) at every level in the dendrogram. The realized deviates, defined as the standardized difference between a fusion value and the mean fusion value, are compared to find significantly large realized deviates. A large realized deviate is significant if the fusion value is at least 1.96 standard deviations from the mean. The upper tail test applies the upper tailed t-test to the fusion values using t-statistic of realized deviate*(n-1) over n-2 degrees of freedom to find the significant increases in fusion value. The bootstrap test performs 1000 trials of randomizing the data, cluster analysis on the randomized data and then compares the random trees to the actual tree. The data matrix was randomized by randomly permuting the order of entries in each column independently. This scrambles the values of each variable for a given cell, thus disturbing any correlations between variables for each cell, but not changing the distribution of values for each variable. The random trees were compared to the actual tree by plotting fusion value vs. number of clusters for the actual data and the distribution of the randomized data with a confidence interval of 1 standard deviation around the mean. Significant departures of the actual fusion values from random, determined by the t-statistic, indicate a number of clusters that are statistically significant. The first level at which significance was found by all three tests is the partition used (indicated in Figures Figures1A,1A, A,2A2A and and3C–E3C–E by bracketing on the bottom of the dendrogram and cut-off linkage distance by a horizontal line).
Finally, in the morphological or physiological clustering we did not detect any obvious bias towards the individual patching the cell, individual reconstructing the cell, or age of the mouse.
K-means clustering was performed for each dataset using Clustan's Focal Point clustering algorithm, with K equal to the number of significant clusters from hierarchical clustering. Focal Point runs the K-means algorithm several times with different initial case orders, which is an important permutation as the K-means algorithm is sensitive to case order. It ranks these solutions based on the minimization of Euclidean distance squared.
Silhouette analysis measures the quality of the clustering by examining the within cluster distances and between cluster distances (Rousseeuw, 1987). The silhouette value of a data point i in cluster A is
where a(i) is the average Euclidean distance between i and all other datapoint in cluster A, b(i) is the smallest average Euclidean distance between i and all datapoints in any other clusters (the average Euclidean distance between i and all datapoints in cluster B, the nearest neighbor to cluster A). The value of s(i) is between −1 and 1. The silhouette width of a clustering is defined as the average of s(i) for all datapoints. Large positive values indicate that the cluster is compact and distinct from other clusters. In addition to calculation of the silhouette width of each clustering, the silhouette width of randomized clustering for each dataset was calculated. The data was randomized as described above for the bootstrap test. The same procedure for cluster analysis was performed on 500 randomized matrices of each dataset and the silhouette width of each clustering was computed. S(random) is the average of the 500 silhouette widths.
To determine which variables are important in distinguishing between the groups identified by cluster analysis, the group mean and standard error was calculated for every variable and the groups were compared for significant differences (Mann Whitney U Test). Additionally, the correlation between the PCs and the original variables was calculated as another measure of which variables are discriminating. Variables highly correlated with the principal components used for cluster analysis are more discriminating. Most of the variables with significant differences between the groups were highly correlated with the PCs (Tables (Tables66 and and77)
To identify subtypes of SOM-positive neurons, morphology and electrophysiology were quantitatively characterized. GFP-positive neurons were patched, recorded from in current clamp mode, filled with biocytin and processed for morphological analysis. 19 electrophysiological variables were measured from current clamp recordings. In addition, 67 morphological variables describing the soma, dendrites, axon and somatic location (relative to the pial surface) were measured from the Neurolucida reconstructions of the patched cells. From a dataset of 59 SOM-positive interneurons from mouse primary somatosensory, visual and frontal cortex, 39 cells had complete reconstructions, 36 cells had high quality electrophysiological recordings, and 16 cells had both complete reconstructions and recordings.
We first explored the statistical structure of the morphological dataset. Cluster analysis using morphological variables was performed with the 39 cells with complete reconstructions, and it revealed five clusters of neurons, grouped into three major branches (Figure (Figure1A,1A, labeled purple, orange and green). One major group encompassed approximately half of the neurons, and could be subdivided into three clusters (purple cells, clusters a, b and c), and was separated by a large Euclidian distance from the other two clusters (orange-d and green-e). In addition, clusters showed no bias with respect to age of the mouse, cortical area, individual patching the cell, or individual reconstructing the cell (not shown).
We further explored the morphological variability of the neurons by plotting their position in principal component space, using the first two principal components as Cartesian axes (Figure (Figure1B).1B). In this analysis one seeks to draw lines that separate each group in a different section of the graph. In this PCA space, the first major group of the cluster analysis (purple cells, corresponding to clusters a, b and c) was clearly separated from the other two branches (orange-cluster d and green-cluster e cells). Within this first group (purple cells), no clear subgroups were observed in PCA space. Clusters d (orange) and e (green) were separated from each other by a second line, with one exception (orange neuron on top). Because of these results in PCA space, although the morphological grouping was significant at the 5 cluster level (black brackets in Figure Figure1A;1A; see Materials and Method and Statistical Validation section below), we chose to group the data into three classes (named groups 1, 2 and 3, and colored purple, orange and green, for the rest of the study).
Examining the morphologies of the cells revealed that group 1 (n=21) corresponded to morphologies traditionally associated with Martinotti cells, a well-characterized subtype of SOM interneurons (Wahle, 1993; Kawaguchi and Kubota, 1996; De Felipe, 2002; Wang et al., 2004). Martinotti cells are found in layers 2/3, 5 and 6 and have a characteristic ascending axon that sends several collaterals to layer 2/3 and layer 1. The axon branches extensively in these layers, particularly in layer 1 where the collaterals branch horizontally, for distances as long as 400μm (see Figure Figure1C,1C, for representative examples; Figure S1A in Supplementary Material, for all neurons)
On the other hand, group 2 (n=11) and 3 (n=7) neurons were quite distinct from Martinotti cells, with very different axonal morphologies. In group 2 and 3 cells, the axon ascended in a very direct course with minimal branching and, in most cases, did not enter layer 1 (Figure (Figure1C1C and Figures S1B,C in Supplementary Material). Axons generally had only one main ascending collateral, or at most, two or three. Initially, we interpreted this sparse axon as a filling or reconstruction artifact, but closer examination of the cells found no evidence of incomplete staining or cut processes. Indeed, endings of their axonal processes did not taper and had apparently normal terminations in the middle of the section (see below). Moreover, group 2 and 3 neurons all had healthy-looking dendrites, without beading, with group 2 cells displaying multipolar morphology and group 3 cells bipolar morphology. Group 2 and 3 cells were found in layers 2/3, 4 and 5.
The combined cluster/PCA analysis indicated that the morphological variance of our database of neurons could be well captured by the hypothesis that they belonged to three major subtypes of cells.
We then sought to use the electrophysiological database to test whether the groups defined by the morphological variables were correlated to electrophysiological groups. For this analysis we used neurons with complete physiological data, regardless of their anatomical reconstruction. Therefore, the two datasets, anatomical and physiological, were not identical but partly overlapped.
When the 36 cells with recordings were clustered by electrophysiology variables, we found significant grouping at the 5 cluster level, by several independent statistical tests (Figure (Figure2A;2A; black brackets; clusters a to e, all neurons with reconstructions colored according to morphological groups described above). The first two clusters (a and b) constituted a group (group 1, purple), in which, remarkably, all of its reconstructed neurons had been previously identified as belonging to the group 1 from the morphological clustering (Figure (Figure1A,1A, purple cells). Besides this first group, the remaining 3 clusters (c–e) contain all neurons previously identified as belonging to the morphological groups 2 and 3 and, moreover, separated these cells in good agreement with these 2 groups. More specifically, physiological cluster c corresponded to morphological group 3 (green) and physiological clusters d and e correspond to morphological group 2 (orange) with one outlier neuron in cluster e that belonged to the first morphological group (purple).
In addition to the cluster analysis, we again explored the statistical structure of the physiological data by plotting the dataset in principal component space (Figure (Figure2B).2B). This PCA analysis revealed a less compact grouping than the morphological dataset. Nevertheless, with two orthogonal lines, one could clearly separate neurons belonging to groups 1, 2 and 3, again, with only one exception, the previously mentioned outlier neuron from cluster e (Figure (Figure2B,2B, purple neuron indicated by arrow). This neuron, interestingly, was close in PCA space to the border between groups 1 and 2.
Inspection of the physiological responses of these three groups of neurons was in agreement with the statistical clustering. Group 1 cells (n=24) displayed a regular spiking firing pattern, non-stuttering and non-fast spiking, with varying degrees of spike frequency adaptation from small to moderate increases in interspike interval (Figure (Figure2C,2C, a and b; average spike frequency adaptation 1.674±0.092; see Ascoli et al., 2008 for definition of electrophysiological variables according to the Petilla nomenclature). For the majority of neurons the action potential (AP) amplitude decayed only slightly over the course of a train (average AP drop 2.619±0.479mV). Group 1 cells also tended to have a lower rheobase (36.042±6.363pA) than group 2 or 3 cells (49.286±10.025pA and 55.000±12.845pA, respectively) and a significantly more depolarized resting membrane potential (Table (Table3).3). On the other hand, group 2 cells (n=7) displayed both regular spiking and stuttering firing patterns. There was a greater spike frequency adaptation (2.244±0.385) and amplitude accommodation than group 1 cells. Group 3 cells (n=5) were regular spiking, with a degree of frequency adaptation similar to group 2 cells (2.657±0.493) and significantly narrower action potentials (Table (Table3,3, AP1 half-width, AP2 half-width). One cell (cell 3, see Figure S2 in Supplementary Material) had an unusual early offset firing pattern of only three to four action potentials per train, and then was silent at a depolarized membrane potential for the remainder of the current pulse. Of the four other group 3 cells, three showed marked amplitude drops in the first few spikes and the other showed amplitude accommodation similar to group 2 cells (Figure S2).
The joint interpretation of the morphological and physiological clustering and PCA analysis was consistent with a simple model of three basic subgroups of SOM neurons. Group 1 is clear: it is composed of Martinotti cells and is clearly delineated in the morphological and physiological clustering dendrograms. Group 2 and group 3 are also distinct in both morphological and physiological clustering, although they are closer to each other in statistical space. Specifically, the close proximity and some overlap between groups 2 and 3 in the principal component space of the morphology dataset versus their complete separation in the principal component space of the electrophysiology dataset suggests that groups 2 and 3 have morphological similarities, although they are physiologically quite different. One could further subdivide the morphological group 1 into three different subgroups (a, b, c; based on the morphological clustering), or the physiological group 3 into 2 subgroups (d and e; based on the electrophysiological clustering). Nevertheless, for simplicity, for the rest of the study we followed the basic three groups (1, 2 and 3), first apparent from the morphological clustering.
These results display a remarkably good agreement between the morphological and physiological classifications, which because they are based on completely independent datasets, confirm each other. In fact, with this basic classification into three groups every neuron is classified into the same groups in both morphological and physiological analysis, with one exception, the outlier purple cell that clusters with group 3 in the physiology, but which is actually borderline in PCA space. Given the tight correlation between morphological group and electrophysiological group for the 16 cells that are in both morphology and electrophysiology datasets, we propose that morphological groups 1, 2 and 3 correspond to electrophysiological groups 1, 2, 3. We do recognize that for a cell belonging to only one of the datasets it is possible that it does not belong to the corresponding group of the other dataset. However, we find this unlikely given that for 15/16 cells the correspondence between morphological and electrophysiological groups was accurate.
We then explored how robust the classification based on cluster analysis actually was, using a variety of statistical analyses. We first used three tests (best cut, upper tail and bootstrap) that are routinely applied to test whether the dataset has a non-uniform structure, an important validation as cluster analysis will always return a clustering tree, regardless of the distribution of the data. Specifically, the best cut test distinguishes the clustering profile of a dataset containing groups from the profile of a dataset with a continuum of points by analyzing the fusion values. Cluster analysis performed on a dataset not containing groups would return arbitrary divisions and thus the fusion values would be continuous, with no large jumps. However, if the dataset does contain groups then at the partitions between the groups large fusion values should be seen. The best cut test finds if there is a linkage cut-off distance at which this fusion value profile occurs. The upper tail test also locates significantly large fusion values. The bootstrap test determines the linkage distance at which the results significantly deviate from random. Indeed, using these three tests, the clusters (as marked by brackets in Figures Figures1A1A and and2A)2A) were found to be statistically significant at 5% level (See Materials and Methods for details). Thus, three independent significance tests all show that the dataset does have a non-uniform structure, and thus could be classified into meaningful clusters.
Following this, we sought to confirm whether the hierarchical classification using Ward's method was reliable, using K-means clustering with K specified as the number of significant clusters from hierarchical clustering. K-means clustering is a top-down clustering method, in contrast to hierarchical clustering which is a bottom-up method. (MacQueen, 1967). While hierarchical classification is inflexible in that cells initially linked stay linked, resulting in a sometimes local optimization rather than global optimization, K-means can reassign cells at a later stage of clustering if it will result in a better optimization. Indeed, the reliability of Ward's method hierarchical clustering for this dataset was proven by the agreement between the two methods (Table (Table4).4). Only a small percentage of cells switched clusters and always stayed within the same group, for example in the 39 morphology cluster cell 27 switched from cluster a to cluster b, both of which are included in Group 1 (Table S3 in Supplementary Material). Such stability between two different algorithms demonstrates that the groups are robust.
In addition, to assess the quality of the hierarchical clusters obtained by Ward's method, we computed the silhouette width of each clustering (Figures (Figures3A,B)3A,B) (Rousseeuw, 1987, See Materials and Methods). Silhouette width is a measure of the appropriateness of each assignment of a datapoint to a cluster. A silhouette width around 0 indicates the datapoint is intermediate between two clusters, silhouette width close to 1 indicates that the point is well assigned and a silhouette width close to –1 indicates that the point is poorly assigned. Constructing a plot of silhouette widths helps reveal the structure of the dataset, and thus to identify natural groups from artificially imposed groups. Using this analysis, we found that the silhouette widths of the clustering by morphology and clustering by electrophysiology confirmed that a natural structure has been found. Such a silhouette width indicates that the distances within clusters are small in comparison to the distances between clusters. Thus, the dataset has a structure containing discrete and compact clusters with clear separation between clusters. In fact, visually inspecting the silhouette plot for each cluster revealed that the majority of cells are well classified.
As a final additional measure of statistical significance, the silhouette width was computed for 500 trials of randomized data, following the same clustering procedure (Table (Table5).5). The average silhouette widths for these randomized trials was drastically lower than for the actual datasets, again demonstrating that the datasets have a strong structure which was disrupted by the randomization.
To determine the quality of classification based on electrophysiological variables, morphological variables and all variables combined, we performed a separate analysis of the 16 cells with both complete reconstructions and recordings (Figures (Figures3C–E).3C–E). This dataset was used to control for differences in sample size and non-overlapping data points between the 36 cell electrophysiology and 39 cell morphology datasets. The clusters found by morphology only, electrophysiology only and both morphology and electrophysiology were in good agreement with the three groups previously identified. Each set of variables had one outlier in the clustering therefore clustering using the combined sets of variables was not an entirely better method. Judging the strength of structure of the three datasets by silhouette width revealed that the clustering by morphological variables has the strongest structure, by electrophysiological variables the next strongest, and by all variables the weakest. The low silhouette width of the combined morphological and electrophysiological variables dataset could be due to the strong correlations within clusters between morphological variables and between electrophysiological variables. With all variables combined, the clusters had a less compact shape due to the multimodality of variables. Thus when clustering by electrophysiology variables only or morphological variables only can discriminate between groups as in our case, we found no particular advantage to clustering by morphological and electrophysiological variables combined.
After establishing the reliability and validity of the classification based on three major groups of SOM interneurons, we explored more carefully the physiological and morphological characteristics of each of them. First, we concentrated on the analysis of the morphological variables that were statistical predictors of each of the three subtypes of neurons (Table (Table3).3). Morphologically, the axon was the most important feature in discriminating between the group 1 and groups 2 and 3, as can be seen by direct inspection of the morphologies. Group 1 cells had typical Martinotti cell axons that extend through layers above the soma and branch extensively, particularly at the pial surface and in the home layer. Group 2 and 3 cells also had an ascending axon, but their axons were much sparser and were contained within a narrow vertical column. Statistically, Group 1 cells had principal component 1 (PC1) values distinct from groups 2 and 3 (Figure (Figure1B)1B) and, indeed, the first principal component of the dataset is highly correlated with 12 axon variables (Table (Table6).6). Six of these variables, axonal node total, total axonal length, total surface area of axon, convex hull axon area, convex hull axon perimeter and convex hull axon surface area, describe the size of the axon in both two and three dimensions. The other variables, i.e. highest order axon segment, axonal Sholl length and node densities and axonal node density describe the degree of axon branching and k-dim axon describes the dimensionality of the axon. Given these correlations PC 1 can be interpreted as a measure of axon size and complexity. Accordingly, there were significant differences for several parameters describing the size, shape and dimensionality of the axon between groups 2 and 3 and group 1 (Table (Table33).
The layers targeted by the axon also differed between group 1 and groups 2 and 3. Group 1 cells extensively branched in layer 1, where Martinotti cells are known to synapse on the apical dendrite tufts of pyramidal cells (Wang et al., 2004; Silberberg and Markram, 2007). In contrast, 13 of the 18 group 2 and 3 cells (72.2%) had axons that avoided layer 1. Some had axons that suddenly terminated at layer 1, while others made an abrupt turn at the layer 1 border and either hooked or continue tangentially along the border (Figure (Figure4A).4A). Interestingly, there was a striking direction preference of this bend at the layer 1 border; in every single case axons turned medially in the slice. Examination of the axons of these cells showed axonal boutons in layer 2/3 (Figure (Figure4B).4B). As mentioned, these axons did not taper and were not cut, terminating in an apparently normal fashion, away from the sectioned surface of the slice (Figure S3B in Supplementary Material). In contrast, the termination of cut axons had larger, more roundly swollen ends and darker ends, and always ended at the very edge of the slice (Figure S3C in Supplementary Material).
The size and shape of the dendritic arbor anatomically distinguished groups 2 and 3, which had overall similar axon morphology. While group 1 was separated from groups 2 and 3 on the PC1 axis, groups 2 and 3 were separated from each other on the PC2 axis (Figure (Figure1B).1B). Examining correlations between the original variables and PC2, revealed that PC2 was highly correlated with several dendrite variables describing the size of the dendrites in two and three dimensions (Table (Table6).6). This possibility is supported by significant differences between the two groups for several variables describing size of the dendritic processes (Table (Table3).3). Group 3 cells had the largest dendritic arbor of the three groups, while group 2 cells had the smallest dendritic arbor, by measures of length and convex hull size. Additionally, group 3 cells had bitufted dendrites, while group 2 cells have multipolar dendrites (Figure S1 in Supplementary Material). Group 1 dendrites were varied, with some multipolar, while others show a preference for extending only upward or downward.
We then focused on the electrophysiological parameters, analyzing both passive and active physiological properties, such as resting membrane state and firing pattern (Table (Table2,2, Figure Figure5).5). There were significant differences among the passive properties of the three groups (Figure (Figure5A).5A). Group 2 and 3 cells were significantly (p<0.025) more hyperpolarized at rest (−70.617±2.385mV and −70.236±2.228mV respectively) than Group 1 cells (−65.336±0.438mV). Group 3 cells had significantly (p<0.025) lower input resistance (256.400±46.427MΩ) than both group 1 and 2 cells (469.583±35.836 MΩ and 475.429±57.726MΩ respectively). In response to hyperpolarizing current, only group 1 cells consistently demonstrated a sag response. This variable was not used for clustering as only 32 of the 36 cells (23 of 24 group 1, 5 of 7 group 2 and 4 of 5 group 3 cells) were given a series of hyperpolarizing current steps. The slope of the regression line fitted to the sag vs. membrane potential was used as the variable for comparison, with sag defined as the difference between the most negative membrane potential during a 1s hyperpolarizing current step and the steady membrane potential in the last 100ms of the step. Group 1 significantly differed (p<0.005 and p<0.01 respectively) from groups 2 and 3 (Table (Table3).3). Some cells showed a depolarized rebound after the hyperpolarizing step, but only one cell in group 1 spiked as a result of the rebound depolarization. As Martinotti cells are known to display a sag response, these results further support that groups 2 and 3 are distinct from Martinotti cells.
We found that the first principal component was highly correlated with the variables measuring the duration, half-width, rise time and fall time of AP1 and AP2 (Table (Table7).7). Therefore it is not surprising that these variables distinguished all three groups (Figure (Figure5B).5B). Specifically, the groups had distinct AP duration, with the values of group 1 (AP1 duration 3.924±0.198ms) being intermediate between groups 3 (AP1 duration 3.000±0.202ms, p<0.05) and 2 (AP1 duration 5.914±0.734ms, p<0.025). In addition, groups 1 and 3 both had faster half-width (AP1 half-width 1.649±0.083ms and p<0.005, 1.356±0.071ms respectively) rise time (AP1 rise time p<0.025, 1.070±0.036ms and p<0.005 0.880±0.037 respectively), and fall time (AP1 fall time p<0.025, 2.853±0.168ms and p<0.005, 2.120±0.166ms respectively) than group 2 (AP1 half width 2.271±0.291ms, AP1 rise time 1.400±0.102ms and AP1 fall time 4.514±0.656ms) while groups 1 and 2 had significantly slower (p<0.025) rise rates (AP1 rise rate 54.265±3.733mV/ms and 48.079±9.906mV/ms respectively) and fall rates (AP1 fall rate 21.812±1.891mV/ms and 17.304±4.749mV/ms respectively) than group 3 (AP1 rise rate 80.337±4.179mV/ms, AP1 fall rate 33.921±2.836mV/ms). The slower rise and fall rate of group 1 was more likely due to group 1 having smaller AP amplitude, while the difference in rise and fall rate between group 2 and 3 resulted from the difference in AP rise and fall time.
In addition, all three groups showed spike frequency adaptation, with group 2 and 3 cells showing a greater degree of adaptation (spike frequency adaptation values group 1: 1.674±0.092, group 2: 2.244±0.385, group 3: 2.657±0.493). Group 2 and 3 cells had greater action potential amplitude as measured from the first two action potentials in a train elicited by twice threshold current injection (AP1 amplitude p<0.05, 61.602±7.865mV and p<0.025, 70.245±2.425mV respectively) than group 1 (AP1 amplitude 55.407±2.548mV) and also had a significantly larger (p<0.025, p<5×10−3) drop in amplitude between these first two action potentials (AP drop 4.822±1.424mV and 6.230±0.911mV respectively) than Group 1 (AP drop 2.619±0.479mV; Figure Figure4A,4A, Table Table33).
Most cells were regular spiking, but three group 2 cells (75, 76, 77) were stuttering, defined as the firing of groups of action potentials that are interrupted by silent periods of varying length (Markram et al., 2004; Ascoli et al., 2008). With only an 8.3% incidence of stuttering cells in the dataset, this is too low a number to conclude that all stuttering cells are group 2 cells, but does indicate that stuttering is a distinctive firing pattern of some group 2 cells. It can additionally be noted that group 2 cells had greater variability in firing patterns, by plotting the dataset in the PC1–PC2 plane (Figure (Figure2B).2B). Groups 1 and 3 were restricted to a much smaller region of the PC1 values, than are Group 2 cells indicating that Group 2 cells have a larger and separate range of values for those variables correlated with PC1 (AP duration, half-width, rise time and fall time).
In this study we used unsupervised cluster analysis to quantitatively classify SOM interneurons based on morphological and electrophysiological characteristics. Using different methods, we find a statistically consistent structure which can be captured well by a simple model of three basic subtypes of cells, already apparent in the first morphological clustering tree (Figure (Figure1A).1A). This does not mean that there are not further subtypes of SOM neurons, but that with the size of our database and methods used in this study, and within the GFP expression pattern of this particular mouse line, at least these three distinct subtypes are apparent. Specifically, the good agreement between clustering based on two independent sets of variables, electrophysiological and morphological, provided strong validation that the groups found by cluster analysis are meaningful and distinct.
Because of the sparseness of the group 2 and 3 cells’ axons we suspected that their axonal morphology was due to incomplete filling or sectioned processes. However, as discussed, we did not find clear evidence of cut axonal terminations in these neurons. While we also did not observe in group 2 and 3 cells any clear signs of incomplete filling, such as dye blown out of the cells or tapering and spotty filling at the end of visible processes, it is still possible that these cells were incompletely filled without these tell tale signs. Another possibility is that these axons were projecting away from the region and we only revealed its local portion. As with other neuronal classes, the projecting axon could be myelinated, and may be difficult to visualize in its entirety. In any case, given that our data was collected from brain slices, we are the first to admit that in vivo evidence is required to conclusively demonstrate the true morphological nature of the axons from groups 2 and 3 cells. Nevertheless, whether their axons are complete or not does not negate our conclusion that they represent different cell types.
Specifically, there are several pieces of evidence in support of group 2 and 3 being real subtypes. First, the axons of group 2 and 3 neurons bend medially, something difficult to reconcile with the possibility that they represented cut group 1 cells, whose axons do not show a characteristic medial bend.
Second, some of the statistical differences in the axonal structures of group 2 and 3 vs. 1 cannot be due to sectioned processes or incomplete filling. For example, differences in the k-dim of axons, or the standard deviation of the tortuosity of axonal nodes, or the axon node densities are not reliant on the size or length of the axon, but instead measure degree of branching and space filling qualities of the axon. If group 2 and 3 cells were the same as group 1, but with cut axons, one would expect these values not to be significantly different across the three groups, because they reflect differences in the topology of the axonal arbor other than size. To test this directly we performed an additional cluster analysis. We excluded all variables related to the size or length of the axon, recalculated principal components and performed hierarchical cluster analysis using these principal components. The following variables were excluded: axonal node total, total axonal length, total surface are of axon, ratio of axonal length to surface area, highest order axon segment, all axonal Sholl variables (6), and all axonal convex hull variables (4). Based on the three statistical tests, the results are significant at the 6 cluster level, and agree well with our results using all morphological variables. There are two outliers, both cells from group 2. One clustered with group 1, and the other with group 3 (Figure S3 in Supplementary Material).
Third, we excluded four cells from our dataset because 3 had cut processes and 1 was not completely filled. These cells resembled Martinotti cells with the ascending portion of the axon cut. Analysis of their axon parameters agreed with this interpretation. The cut Martinotti cells had a mean total axon length significantly larger than groups 2 and 3, but slightly smaller than group 1 (p<0.005 for group 2 and p<0.01 for group 3). These cut cells were also significantly different from group 2 (p<0.005) and group 3 (p<0.01) with respect to axonal node total, total axonal length, total surface area of axon, and was not significantly different from group 1 with respect to the same variables. The cut Martinotti cells also had the same axon branching properties as group 1 cells, highly branched and space filling, distinct from the sparse axon branching of groups 2 and 3 (again, not significantly different from group 1 with respect to axonal node total, k-dim axon, stdev of tortuosity of axonal node, axon node density, p<0.025 for groups 2 and 3 with respect to those same variables).
Finally, in some slices more than one GFP cell was patched, filled and reconstructed. In most cases, all cells were Martinotti neurons. However in two slices there was both a Martinotti cell and a group 2 cell (cell 28 and 29 together and cell 52 and 53 together). Since the axon of the Martinotti cell was large and complete, and the cells were located in a similar section of the slice, it is unlikely that only the axon of the group 2 cell was selectively cut. Moreover, cells 52 and 53 even had overlapping axonal arbors in some portions and there was not evidence of incomplete filling.
Further evidence that group 2 and 3 are real comes from their electrophysiological properties. Specifically, clustering by electrophysiology agreed with clustering by morphology as both separated the Martinotti, group 2 and group 3 cells, providing an independent control for the possibility that the group 2 and 3 cells were poorly filled or reconstructed group 1 cells. The AP durations of the cells are long but consistent with literature values for recordings from juvenile mouse slices at room temperature. We attribute these long APs observed to two reasons. Firstly, juvenile animals have slower action potentials than adults and secondly APs are slower at room temperature than at 35–37°C (Ali et al. 2007). SOM+cells in P13–P16 rat slices recorded at room temperature had AP1 duration of 3.69±0.53ms and AP2 duration of 4.32±0.55ms, similar to our Martinotti cells (Wang et al., 2004). The AP half-width of our cells is also in line with room temperature recordings in juvenile rat slices from bitufted interneurons, some identified as Martinotti cells (Ali, 2007), and from regular spiking interneurons (Cauli et al., 1997). AP durations reported for groups 2 and 3 are not unreasonable given the ranges in the literature for recordings in slice from juvenile animals at room temperature. The RMP and input resistance are in a healthy and normal range for all cells. Finally, as mentioned, group 1 and groups 2 and 3 differed statistically in their sag.
The fact that Martinotti cells comprise a subtype of SOM interneurons was already known, thus the separation of the Martinotti cells is an external validation of the results. Finally, the distinction between group 2 and 3 was preserved in both morphological and electrophysiological clustering, suggesting that the non-Martinotti cells are comprised of at least two subtypes. These two groups do not appear to be any of the SOM subtypes previously identified in mouse cortex (Wang et al., 2004; Halabisky et al., 2006; Ma et al., 2006), although they do have some similarities to neurons from the X94 somatostatin positive line, such as avoidance of layer 1 (groups 2 and 3) and stuttering properties (group 2). However, X94 cells have highly branched axons, contained mostly in layer 4. Thus, X94 cells and our group 2 and 3 cells target different layers. Additionally X94 cells are quasi-fast spiking, with much shorter ISI than group 2 or 3 cells. For the four electrophysiological subtypes identified by Halabisky et al., two of the subtypes described have a similar level of spike frequency adaptation to our cells. In addition, like group 3 cells, all four subtypes they identified have a low input resistance with mean input resistance for the four subtypes between 205 and 224MΩ. However, in this study the sEPSC kinetics in addition to passive and active membrane properties were important in distinguishing between the four subtypes and morphologies of the cells were not recovered. Therefore, comparison of our results to these subtypes is difficult without sEPSC or morphological comparison. In futures studies it will be interesting to determine the molecular phenotype of these SOM- positive interneurons and identify any differentially expressed calcium binding proteins or neuropeptides.
The use of juvenile mice (P10–P18) raises the question of whether the subtypes identified in this study are representative of cell types in adult mice. A comparison of the properties of layer 4 neocortical interneurons in slice from juvenile rats (P18–22) and adults rats (Ali et al., 2007) indicates that subtypes identified in juveniles likely are representative of adult cell types. This study found that there was no significant difference in axonal arbor size between the juvenile and adult. Correlations between firing properties and synaptic properties were also conserved from juvenile to adult. The only significant difference identified was that in adult animals, action potential, EPSP and IPSP half-widths were narrower and AP-AHPs were deeper than in juvenile animals. Given that the relationship between firing properties and synaptic properties is not dependent on age, findings in juvenile animals are relevant. In particular the identity of different classes of interneurons described in the study remains the same between juveniles and adults. For example the bitufted interneurons had broad action potentials relative to other interneuron types in both juvenile and adult animals. As all cell types had slower events in juvenile tissue, the differences between groups remain unaffected by age.
Given the increasingly sophisticated understanding of SOM interneuron development, it is imperative to discuss whether there is any evidence for different subtypes of SOM interneurons based on this knowledge. Of particular relevance is the clear relation between the neuronal cell type and spatio-temporal anlage of fate determination (Wonders and Anderson, 2006; Dasen and Jessell, 2009). SOM interneurons, along with the majority of neocortical interneurons, originate from the medial ganglion eminence (MGE) (Xu et al., 2004; Wonders and Anderson, 2006), although it has been so far difficult to pinpoint exactly the precise spatio-temporal profile of their precursors. By using inducible genetic fate mapping, recent studies have addressed whether different interneuron types originate at different points in development. Based on molecular content and electrophysiology, distinct groups have been reported (Miyoshi et al., 2007; Sousa et al., 2009). By fate mapping Olig-2 expressing precursors (restricted to MGE where Olig-2 expression is high) at different times in development, it was found that in mouse somatosensory cortex barrel field SOM-positive, calretin (CR) negative interneurons make up 30% of early populations (precursors labeled at E9.5–10.5), but their generation drops as development continues, comprising almost none of the late population (precursors labeled at E15.5) (Miyoshi et al., 2007). In an analogous study fate mapping Nkx2.6 precursors, a similar temporal pattern was observed for SOM-positive CR-negative interneurons. In addition they found that SOM-positive CR-positive interneurons make up the majority of the late population (precursors labeled after E12.5; Sousa et al., 2009). The SOM positive CR positive interneurons have multipolar dendrite morphology and regular spiking patterns, while the SOM positive CR negative interneurons are identified as Martinotti cells. Therefore it is possible that the three groups identified in this study originated at different developmental points, or perhaps in different spatial positions. We would therefore venture the hypothesis that the lack of consensus in identifying the origin of SOM interneurons reflects their basic heterogeneity, so their classification into three (or potentially more) subgroups could help better delineate the spatio-temporal origin of each subgroup.
The very distinct axon morphology of group 2 and 3 cells may be the result of a common developmental response to molecular cues, and is thus congruous with the groups originating at different times or from different locations in the MGE. The avoidance of layer 1 by group 2 and 3 cells, while group 1 cells have extensive branching in layer 1 could be due to either opposite responses to the layer 1 molecular cues or attraction of group 2 and 3 cells to axon guidance cues in layer 2/3. Such a feature can be seen on some of the layer 4 short axon cells drawn by Lorente de No (1922).
The distinct morphology and physiology of each of the three types implies that each could have a different functional role within the circuit. First, the differently shaped dendritic arborizations mean that the groups have different laminar inputs. In addition, differing axon morphologies imply that downstream targets of the groups are distinct. In particular, Group 2 and 3 cells’ axons avoid layer 1 while Martinotti cells branch extensively in this layer. The modulation of dendritic excitability in layer 1 is important as connections from higher cortical areas are received in layer 1 (Coogan and Burkhalter, 1993). Through a connected network of pyramidal cells and Martinotti cells, excitation of a single pyramidal neuron can recruit Martinotti cells which in turn cause inhibition of the apical tuft dendrites of neighboring pyramidal cells (Kozloski et al., 2001; Kapfer et al., 2007; Silberberg and Markram, 2007). This pathway is extremely sensitive, with a supralinear increase in recruited Martinotti cells as the number of pyramidal cells excited increases (Kapfer et al., 2007).
The medial directional preference of group 2 and 3 axons suggests that this asymmetry could have a functional purpose. This asymmetry resembles other incidences of neurons with asymmetric processes with potential functions. For example, retinal OFF direction selective ganglion cells have asymmetric dendrites with a dorsal to ventral direction preference. This asymmetry translates to an asymmetric and selective receptive field which allows these cells to detect only upwards movement (Kim et al., 2008). These retinal ganglion cells have a distinct molecular identity and “tile” the retina, such that each distinct ganglion cell type covers the retina with minimal overlap (De Vries and Baylor, 1997).
In the neocortex there are also examples of asymmetric dendrites and axons. In macaque monkeys, Meynert cells in layer 6 of primary visual cortex have asymmetric basal dendrites, extending up to 600 μm and aligned parallel to the ocular dominance column (Winfield 1983; Livingstone, 1998). These cells have a direction selective receptive field based in part on the dendritic conduction delays due to asymmetry of the basal dendrites (Livingstone, 1998). In primary somatosensory barrel neurons in layer 4 have asymmetric dendrites and axons that are almost always directed toward the center of the barrel (Woolsey 1975; Egger et al., 2008). It is postulated that this allows the dendrites to receive maximal input from thalamocortical afferents, while the localization of the axon provides minicolumns of recurrent excitation within the barrel. Distinct minicolumns could each process a receptive field for a different quality of whisker movement. A few barrel neurons have axons directed to the neighboring barrel, providing interbarrel connections and allowing for multiwhisker receptive fields (Egger et al., 2008). Such examples of structure so clearly relating to a particular function could hold true for group 2 and 3 cells whose axons appear to have a very selective target, given the sparsity and directional preference of the axon.
The checkered history of neuronal classification based on qualitative descriptors has led to a multitude of classifications. Unsupervised classification methods appear ideal to objectively identify the subtypes of neurons present in biological circuits using quantitative criteria obtained by blind algorithms. At the same time, it is possible that supervised methods, using available information in the dataset, could outperform clustering algorithms for certain neuronal classification tasks (L. Guerra, L. M. McGarry, P. Larrañaga, V. Robles, C. Bielza and R. Yuste, in prep.), so further methodological research appears necessary.
As cluster analysis becomes a more widely used tool for defining cell types, it is important to note that it should used in such a way as to optimize its power of objective and quantitative classification, but while also controlling for its weaknesses. One problem of cluster analysis is that it always returns a classification of the input data. Thus it is important that the results of cluster analysis are thoroughly validated for statistical significance and evidence of natural groups, rather than artificial or arbitrary divisions. For example, the quality of the clusters should also be examined with silhouette analysis as well as agreement between classifications based on multiple different independent sets of parameters. Additionally the weaknesses of an algorithm being used can bias the results. Thus, comparison to an alternate clustering algorithm, preferably a complementary algorithm such as bottom up hierarchical clustering paired with top down K-means clustering, tests for robustness of the classification. All of these approaches were used in this study to assess our cluster analysis results, allowing us to derive more convincing evidence for our conclusions. Finally, in this study with Wards’ hierarchical clustering, we have encountered no particular advantage to use a multimodal cluster analysis (i.e., one that included physiological and morphological features) over using just one (for example, the morphology, Figure Figure3).3). While this could reflect the structure of our data, it should be kept in mind as a warning note for future studies using Wards method of multimodal datasets.
Using cluster analysis raises a fundamental question: what should be considered a distinct cell type? As cluster analysis will generate continual division, it becomes necessary to define an endpoint that is both meaningful and useful in addition to being statistically significant. Part of the problem is the uncertainty about the number of neocortical interneuron types, with estimates of as many as 100 cell types (Lorente de No, 1922; Stevens, 1998). There is still much controversy with regards to definition of a cell type, but several recommendations made by recent reviews have begun to address the issue, proposing some objective criteria, such as tiling of surface (Mott and Dingledine, 2003; Masland, 2004; Migliore and Shepherd, 2005) or their transcription factor expression (Yuste, 2005). The eventual goal is to define functionally distinct groups, most often identified through a combination of characteristic morphology, electrophysiology or molecular profile. The continuation of quantitative classification of interneuron subtypes, coupled with new findings on both the developmental origin and functional role of these subtypes will make this goal come closer to realization. Such a complete description of interneuron subtypes is necessary to achieve the final goal of mapping out the cortical circuitry.
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
The Supplementary Material for this article can be found online at http://www.frontiersin.org/neuralcircuits/paper/10.3389/fncir.2010.00012/
We thank Luis Guerra for advice, members of the laboratory for reconstructions, help and comments and anonymous reviewers for helpful suggestions. Supported by the Kavli Institute for Brain Science and the NEI (EY11787).