article: Every neuron is part of a network, exerting its function by transforming multiple spatiotemporal synaptic input patterns into a single spiking output. This function is specified by the particular shape and passive electrical properties of the neuronal membrane, and the composition and spatial distribution of ion channels across its processes. For a variety of physiological or pathological reasons, the intrinsic input/output function may change during a neuron’s lifetime. This process results in high variability in the peak specific conductance of ion channels in individual neurons. The mechanisms responsible for this variability are not well understood, although there are clear indications from experiments and modeling that degeneracy and correlation among multiple channels may be involved. Here, we studied this issue in biophysical models of hippocampal CA1 pyramidal neurons and interneurons. Using a unified data-driven simulation workflow and starting from a set of experimental recordings and morphological reconstructions obtained from rats, we built and analyzed several ensembles of morphologically and biophysically accurate single cell models with intrinsic electrophysiological properties consistent with experimental findings. The results suggest that the set of conductances expressed in any given hippocampal neuron may be considered as belonging to two groups: one subset is responsible for the major characteristics of the firing behavior in each population and the other is responsible for a robust degeneracy. Analysis of the model neurons suggests several experimentally testable predictions related to the combination and relative proportion of the different conductances that should be expressed on the membrane of different types of neurons for them to fulfill their role in the hippocampus circuitry.
  Any given neuron in the brain is part of a network, in which it exerts its action by transforming the input it receives into an output. This function is specified by the particular shape and passive electrical properties of the neuronal membrane, the composition and spatial distribution of ion channels across its processes, and the functional properties of the synaptic inputs themselves. During development and during the entire lifetime of a neuron, its input/output function is adapted to realize ongoing refinement of the function of the neuron and circuit, or maintain functional robustness in the face of constant protein turnover or an evolving pathological condition. Such adaptability of individual neurons can be achieved through a myriad of dynamic mechanisms, including structural, intrinsic, and synaptic plasticity. A direct experimental evidence for these mechanisms is the high variability observed for the current generated by specific types of ion channels measured across individual neurons, from either a homogeneous population or different cell populations (e. g. [1]). The mechanisms responsible for this variability are not well understood, although there are clear experimental and modeling indications that correlation and degeneracy among a variety of conductances can be involved [2,3]. The phenomenon of degeneracy allows the possibility, for a complex biological system, to perform the same function using structurally different elements [4]. In the context considered in this paper, it refers to the robust and tunable adjustment of a neuron’s firing properties [5]. For example, a neuron can be tuned to perform a given function by expressing in the membrane a specific set of conductances with a specific dendritic distribution (Migliore (2003) ); degeneracy can result in this tuning being robust, by implementing the same function with many different configurations of the same set of conductances. This property has been systematically studied in crab stomatogastric ganglion neurons [2,6] and in Globus Pallidus neurons of the rat [7]. In the present study, we investigate this issue for neurons of the hippocampal CA1 region. These neurons are important because they have a critical position as the main output stage of the hippocampal circuitry [8]. The hippocampal CA1 pyramidal neurons, in particular, exhibit a peculiar ensemble and distribution of conductances (reviewed in [9]), subject to significant changes following activity-dependent biochemical processes, such as activation of protein kinase A and C, or Ca/calmodulin dependent kinase II [10,11,12], pathological conditions (e. g. [13,14]), or traumatic brain injuries [15,16]. There must then be an extremely robust compensatory mechanism in these neurons, or in the network, which maintains or re-establishes the physiological activity within an operation range, in spite of a potentially large change in its intrinsic properties or synaptic input. Here we study the mechanisms of robustness of intrinsic properties by using a unified data-driven workflow and open source analysis and simulation tools. From a set of experimental recordings and morphological reconstructions, we implemented many morphologically and biophysically accurate models for CA1 pyramidal neurons and interneurons, with intrinsic electrophysiological properties constrained by and consistent with the experimental findings. The results indicate that a few currents need to be expressed at a relatively stable level, whereas others can be expressed within a much wider range. The analysis of the model neurons suggests many specific experimentally testable predictions on the combination and relative proportion of the different ionic conductances, and their relationship to robustness of intrinsic properties.
  To implement a set of data-driven neuron models, we start from a set of morphological reconstructions of neurons and somatic voltage traces obtained from in vitro slice preparations of rat hippocampal tissue to use as constraints (see Methods). In Fig 1 we show several examples of the 34 morphologies used in this work (19 pyramidal cells and 15 interneurons), superimposed on a rat hippocampal slice stained for parvalbumin for illustrative purposes. A total number of 1456 experimentally obtained somatic voltage traces for a range of stimulation protocols were used in the optimization pipeline to constrain the models (see Methods). Collections of traces for individual neurons were manually assigned to four electrical types (e-type), according to the firing pattern exhibited during increasing somatic current injections [18], and using the classification proposed in the Petilla convention [19]. The 832 traces from pyramidal neurons, with an increasing inter-spike-interval (ISI), were all classified as continuous accommodating cells (cAC). For interneurons, 240 traces were classified as cAC, 160 traces as bursting accommodating cells (bAC), and 224 traces, whose firing rate is constant, as continuous non-accommodating cells (cNAC). Typical examples illustrating the physiological variability observed for these e-types are shown in Fig 2. A more quantitative analysis and comparison of their features will be presented elsewhere (Bologna et al. , manuscript in preparation). Different pyramidal neurons (Fig 2, pyr cAC) exhibited significantly different responses to the same input. For example, a near-threshold 0. 4 nA somatic current injection may or may not generate a few action potentials, whereas a 0. 8nA input can result in a 2-fold range for the number of elicited action potentials (APs) (Fig 2, pyr cAC, blue traces). Interneurons classified as cAC also exhibited a large inter-cell variability, with different cells responding to the same stimulus with a wide range of spike patterns, such as tonic firing (Fig 2, int cAC plots, cell 970428A1), stuttering (cell 970509HP2), and depolarization block (cell 980205FHP). The other two interneuron e-types, bAC and cNAC, also exhibited a large variability among different cells (Fig 2, bottom plots). This variability can be the result of different morphologies and/or a different density and distribution of the conductances expressed on the membrane of the different neurons. In the following sections, we will explore in more detail this issue by implementing and analyzing cellular level models that are able to reproduce these results. For each e-type (see S1–S4 Tables and Methods), a set of electrophysiological features were extracted from all voltage traces belonging to that e-type. All the pyramidal cell morphologies were used to implement cAC models, whereas interneuron morphologies were used to obtain cAC, cNAC, and bAC models following the known firing behavior of each type of morphology (see legend of Fig 1 and S5 Table). Features and morphologies were then used to obtain a set of optimized models for each e-type, using a heuristic parameter optimization process that employed multi-objective genetic algorithms. Each optimization run (see Methods for details) returned a number of viable “individuals”, each one with a specific ensemble of peak ion channel conductance and passive properties consistent with the chosen “objectives” (i. e. a set of experimental features). As a cost function for the optimization process we used a score defined by the total error associated with each individual, calculated as the sum of the absolute deviations of model features from the experimental mean, in units of the experimental standard deviation (sd) obtained for the value of each objective. A score = 0 would correspond to an individual with all parameters equal to the average value of the corresponding experimental electrophysiological feature. The total error thus gave an idea of how good the individual was in representing the neuron’s overall expected behavior under a series of 400 ms long somatic current injection steps. The final choice to accept an individual as a plausible representation of a given e-type was based on the error obtained for each objective. An individual with a sd<2 for all objectives was considered acceptable. Typical optimization results for pyramidal and interneuron cAC e-types are shown in Fig 3. Traces obtained for different somatic current injections from three individuals (Fig 3, traces on top left graph of each panel), showed that the optimization process was able to take into account the experimental variability. Different individuals exhibited significantly different responses to the same stimulus, as in the experiments. The evolution of the total score as a function of the number of generations in the optimization process (bottom graph in each panel), showed that the optimization converged nearly monotonically in relatively few iterations, having reached a relatively stable minimum within approximately 60 generations. The list of objective scores for the best individual in each case (Fig 3 right graph in each panel) showed that for most features (n = 60 for pyramidal cells and n = 47 for cAC interneurons, see S1–S4 Tables) the associated error was below 2 sd. Similar results were obtained for the optimizations of bAC and cNAC interneurons (see individual optimization files at https: //collab. humanbrainproject. eu/#/collab/18565). Taken together, these results show that the overall optimization process is a robust way to obtain a number of biophysically accurate neuron models of hippocampal CA1 pyramidal cells and interneurons, which are able to reproduce many of the properties observed experimentally in different types of neurons. A more direct comparison between experimental and modeling traces for the different e-types is shown in Fig 4A, revealing a very good qualitative agreement between the modeling results and experimental traces. The optimization enabled the production of models that correctly reproduced many characteristics of the firing patterns, such as the strong accommodation observed in cAC interneurons (Fig 4A, cAC int @0. 4nA), the high firing frequency of bAC interneurons at the beginning of a current injection (Fig 4A, bAC @0. 6nA), and the progressive reduction in the AP amplitude during the first part of stronger stimuli (Fig 4A, bAC @1nA). The pyramidal cell models also exhibited a typical property often observed experimentally in this type of cells, i. e. the decrease in the peak amplitude of an AP backpropagating into the apical dendrites [20]. This effect has been shown to depend on the high density of A-type potassium channel in the apical dendrites [21], but not all CA1 pyramidal neurons exhibit this effect [22,23]. It is important to note that this feature was not used to constrain the optimization but, interestingly, the optimized models were able to reproduce it, as shown in Fig 4B, for a few cases using morphologies from both young adult (cells 050921AM2, and 990803) and P14-23 animals. The dichotomy in AP backpropagation observed in the experiments [22] was also reproduced by the model neurons, with the AP amplitude either strongly decreasing beyond ~150 μm from the soma or limited to ~50% of the maximum, with very few cases in between. Taken together, this comparison between experiments and models at the individual trace level, suggests that the optimization process was able to correctly capture and explain both intra- and inter-cell variability in firing behavior in terms of different combinations of active and passive membrane properties. An indication of how the optimized models may capture the variety of experimental input/output properties can be drawn from Fig 5, where the number of spikes for each e-type was plotted against the somatic current injection, for experimental (blue lines) and modeling traces (red lines). In all cases, experimental traces exhibited a rather large inter-cell variability in the number of spikes elicited by any given input current. It is quite common to see up to a ~5-fold difference in the number of spikes elicited in different cells under the same current injection. In most cases, the models were in quantitative agreement with the average number of spikes generated as a function of the input current (Fig 5, insets, Mann Whitney Rank Sum test p>0. 05 in all cases except for 1nA injection in pyramidal neurons). With the set of data-driven neuron models obtained for each e-type, we can now analyze how different combinations of peak conductances can result in models able to reproduce equally well the firing properties observed experimentally under different current injection steps. The optimization process generates many of these models (termed “individuals”) because of ion channel degeneracy [5]. As discussed in the Introduction, this phenomenon is thought to allow a neuron to adjust its firing properties in a robust and tunable manner. To obtain further insight into on how degeneracy is achieved in hippocampal CA1 neurons, we analyzed all the individuals obtained from the optimization runs. For each optimization run, the 10 best individuals were considered based on their total score (see Methods). Note that these individuals were obtained from the same morphology with different channel densities. In Fig 6, the value of the optimized parameters, normalized to the maximum value chosen for each conductance, were plotted for each optimization run (10 individuals for each run, opt id). For clarity, in each graph the values obtained for any given parameter were placed on the Y-axis according to the corresponding average value calculated from all optimizations. In this way, the bottom rows in each graph correspond to parameters with an average low value whereas top rows correspond to parameters with higher values. Furthermore, parameters that were relatively stable across all optimizations (i. e. with a sd<0. 2) for any given e-type are highlighted using a red label in the y axis. For pyramidal cells (Fig 6, pyr cAC) the most stable parameters were some of the passive properties, Ih, KM, Calcium, and Ca-dependent K currents. Interestingly, we noted that whereas passive properties were consistently optimized with a stable value across the optimizations for all e-types (Fig 6, see top rows in all graphs), conductances were shown to be somewhat different depending on e-types. For example, for interneurons, Ih, somatic KM and dendritic KDR were the most stable for all e-types, whereas dendritic KA was stable for cAC and Cagk for cNAC. These results suggested that each e-type has specific active properties that may be particularly important to obtain the appropriate firing pattern in response to a given input. While these properties need to be well constrained for each e-type, degeneracy can be achieved by combining the other conductances in a relatively large number of ways. The functional consequences of this situation will be discussed below. To explore whether a cell’s morphology can also be related to degeneracy, we fixed the peak conductance values to those found for the best overall individual (obtained for morphology oh140521_B0_Rat_idA) and calculated the total error by using different morphologies. The results are shown in Fig 7A. We found that the total error using the same set of conductances on different morphologies was within the range obtained for each cell’s optimization for 10 out of 16 morphologies. For these cases, there was no correlation between the total error and the main morphological properties, such as soma area, total cell volume, or number of sections (Fig 7B, Spearman correlation, p>0. 05 in all cases). These results suggest that degeneracy can also be obtained using different morphologies equipped with identical peak channels conductance. A deeper analysis of this issue however was not further considered in this work. For a more detailed analysis of the configuration of peak conductance values for all models, we first considered the results for pyramidal neurons. In Fig 8A we show a typical distribution of normalized values obtained for membrane properties where optimizations yielded a relatively narrow range (somatic KM, Ih, and Ra), or a wider range of values across individuals (dendritic Na). Note that two of the conductances with a narrow distribution are, in pyramidal CA1 neurons, the dominant factors in controlling major properties such as excitability and accommodation (KM, reviewed in [24]), and synaptic integration (Ih, [25]). The paramount importance of these two types of conductance for reproducing the experimental traces, suggested by their value lying in a narrow range across individuals, emerged from the optimization process without any specific constraint. An insight on degeneracy in these neurons can be obtained by considering correlation between parameter pairs. In most cases, we found no statistically significant correlation (see S6 Table for the Spearman correlation coefficients). However, for several cases a significant correlation between selected parameters was found (S6 Table, grey cells). The conductance which was most correlated with others was Cagk, a Ca- and voltage-dependent K+ conductance that is one of the major determinants for accommodation in these neurons. The inverse correlation with the KM is particularly interesting, since it supports the experimental finding that these channels operate in combination to control intrinsic hyperexcitability [26], and modeling results suggesting how they must both be involved to obtain a strong accommodation [27,28]. To explore the configuration of the conductances in a more qualitative and intuitive way, we arranged a radar plot of the conductances most correlated with Cagk (Fig 8B), and one of those showing little variability (in this case the reversal potential of the leakage current in the dendrites, e_pas d). The different individuals were sorted with respect to Cagk (Fig 8B, thick black line) and, for clarity, we plotted only 40 of the 160 individuals. The highly jagged and intermixed lines represent the different peak conductance type and value for different individuals giving equally good representations of 60 electrophysiological features experimentally observed in these neurons (see S1 Table). Examples of model traces from a few individuals (all obtained with a 0. 4nA somatic current injection) displayed the same number of spikes obtained with very different channel configurations. The number of spikes elicited for each individual is plotted in Fig 8C. These results give a clear indication that degeneracy in CA1 pyramidal cells can easily emerge from many different combinations of many, but not all, channels. The reason for the lack of pairwise correlation between most parameters does not exclude that the parameter space may be shaped by higher order correlations that can be ultimately responsible for degeneracy. However, a full quantitative study of higher order correlations was outside the scope of this study. The results obtained for interneurons are shown in Fig 9. In this case, to allow an easier comparison of the parameters among the different e-types, individuals were sorted according to the somatic Na conductance (Fig 9, thick black lines), which was among the most correlated with all the others (see S7–S9 Tables). The models suggest a few distinct differences among the different e-types. Note, for example, the distribution of values obtained for the peak conductance of dendritic KDR or KA in the various e-types (Fig 9, dark red and blue lines, respectively), or the difference in the overall values of dendritic Na (Fig 9, orange lines) between cAC and cNAC. In general, however, the distribution of values were analogous to those obtained for pyramidal cells, with each individual characterized by a highly variable combination of values for many conductances. Finally, one important factor in determining the firing characteristics of different neurons, in addition to a substantial change in morphology [29] and/or gene expression profile [30], is the relative proportion with which specific channels are expressed on the membrane. For this reason, from the optimized models we calculated the relative contribution of each channel in each e-type, by considering the average value of each peak conductance calculated across all individuals. The results are presented in Fig 10A. In all cases, we found that Na, KA and KDR could account for most of the channels expressed on the membrane. Interestingly, each e-type showed a distinct proportion of these channels, with axonal Na channels playing a relatively large role in all e-types, axonal KA being relatively more important in pyramidal neurons than in interneurons, and dendritic KDR being significantly higher in cNAC e-types. An analysis of the relative level of each conductance in the various e-types (Fig 10B) also showed significant differences in several cases (Pairwise Multiple Comparison Procedure, p<0. 05). From the results it is clear, for example, that dendritic Na should be higher in pyramidal cells than in any type of interneuron, cAC interneurons should have a higher dendritic Na among interneurons (Fig 10B, dark blue squares for Na d), and that the axonal KM is essentially independent from cell type. In summary, these results suggest the experimentally testable prediction that different e-types can be characterized by a different combination of the same set of conductances.
  It has been shown that any individual neuron can express a distinct combination of many channel types [30] determining its electrical properties [31]. Furthermore, several seminal papers demonstrated that each cell type could exhibit specific correlation between channels expression [32], which may emerge from a homeostatic rule [2]. The overall picture is one in which many different conductances coincide to produce the electrophysiological patterns that characterize the operating range of any given population of neurons, and they do so in such a way to compensate for relatively large changes in individual channel density or synaptic connectivity [33]. The robustness of this mechanism relies on degeneracy [4], which can be practically implemented through a large and flat parameter space for channel conductance. This issue has been studied in the crab pyloric neurons [3], stomatogastric ganglion neurons (e. g. [2,6]), in the Globus Pallidus neurons of the rat [7]. The presence of degeneracy had yet to be studied in hippocampal neurons. Two recent modeling studies, in the mouse corticospinal neurons and motor cortex, have explicitly shown how degeneracy in cortical neurons can work to implement some electrophysiological features but not others [34], and that degeneracy can also generate multitarget routes from pathological to physiological network dynamics [35]. The first finding was particularly relevant for our study, and it was among the reasons why we choose not to include the voltage between spikes among the optimized features. Its accurate reproduction would have required us to additionally optimize channel kinetics, which was not within the scope of this work. The analysis of the modeling results presented in this paper provides many experimentally testable predictions on the possible co-regulation of ion currents in hippocampal CA1 neurons. Correlation between pairs of specific conductances has been found for cells in the stomatogastric ganglion of the crab (STG, [32]) and in the pyloric network of the spiny lobster [36]. These experiments found that several pairwise correlations between the same conductances can be present in different type of cells, but no cell type showed conductances with the same set of pairwise correlations. Our optimized models confirmed this result also for the hippocampal CA1 neurons. The models also confirmed pairwise correlations already observed in STG, such as that between KA and Ih, Na, KDR, and Cagk, and between Na and Cagk. Like in the STG, these correlations were observed in different combinations among different cell types. It is important to stress that the optimization process did not bias the parameter values against each other. Correlations thus emerged naturally from the optimization process, and reflected a better reproduction of the experimental features. The models predict several additional pairwise correlations between conductances (see S6–S9 Tables), which are specific for each e-type. All predictions can be tested experimentally, by directly measuring and comparing peak ion currents or (better) channel densities in different neurons or by a genetic perturbation of channel expression [36,37]. A limitation of this work is that the optimization process was not able to generate a population of models reproducing the very large experimental variability. The reason for this effect is that, in this work, we choose to optimize the different e-types using for each feature the average and standard deviation calculated from all traces, rather than independently optimizing models constrained by traces from an individual cell. A partial explanation for this choice was the limited availability of experimental data on individual cells. Nevertheless, we think that these results offer a significant improvement on the current state of the art, and a necessary step towards building a full-scale cellular model of the rat hippocampus CA1 circuit (Romani et al. , in preparation). Another experimentally testable prediction of the models is that each type of cell should have a small number of channel types that would be expressed at the same density in the same neuronal population. There is already some experimental indication that this is the case for STG cells in the crab [1], where it has been found that KDR is relatively constant among the lateral pyloric neurons of different animals, whereas KA and Cagk varied more than threefold. In this study, we found that passive properties, KM, and Ih were among the most stable intrinsic membrane properties in any given neuron population, together with dendritic KDR for interneurons. The models also predict that a different combination of axosomatic Na, KA, and KDR channels may dominate the distribution of channels on the membrane of a neuron belonging to a given e-type. This is also experimentally testable, by directly measuring the density of the different channels expressed on the membrane of different type of neurons. Our analysis suggests a physiological plausible explanation for why single channel mutations can have more or less pathological consequences. A clear example stands out for KM and Ih channels in pyramidal cells. We found that these channels must be expressed with a relatively stable density; they do not appear to contribute to degeneracy. This may explain why specific mutations of KM channels can result in neonatal epilepsies in humans [38], or why the decrease in Ih caused by experimental models for temporal lobe epilepsy can result in major changes in the electrophysiological mechanisms related to cognitive functions [39]. Finally, the modeling effort presented and discussed in this work is part of a larger modeling workflow currently underway in the framework of the EU Human Brain Project (https: //www. humanbrainproject. eu/en/), with the main goal to implement a cellular data-driven model of the entire hippocampus. The Hippocampus is a complex brain structure, deeply embedded into the temporal lobe, with a paramount importance for higher brain functions such as learning and memory, and spatial navigation, and is involved in several major brain diseases. In spite of intensive experimental and computational studies, the mechanisms underlying these functions (and dysfunctions) are still poorly understood. A model implementation and analysis at the cellular level may pave the way for a deeper understanding of the diverse and complex functions of this brain region, and of its levels of organization. One of the major steps towards this goal is the implementation of morphologically and biophysically accurate single cell models for the main neuronal populations, equipped with a set of axonal, somatic, and dendritic currents consistent with many experimentally measured electrophysiological features, in such a way as to be able to capture the main I/O properties observed experimentally. Here we have used a general, robust, and flexible tool able to produce, using reasonable computational resources, ensembles of this type of models for CA1 pyramidal cells and interneurons.
summary: The peak conductance of many ion channel types measured in any given animal is highly variable across neurons, both within and between neuronal populations. The current view is that this occurs because a neuron needs to adapt its intrinsic electrophysiological properties either to maintain the same operative range in the presence of abnormal inputs or to compensate for the effects of pathological conditions. Limited experimental and modeling evidence suggests this might be implemented via the correlation and/or degeneracy in the function of multiple types of conductances. To study this mechanism in hippocampal CA1 neurons and interneurons, we systematically generated a set of morphologically and biophysically accurate models. We then analyzed the ensembles of peak conductance obtained for each model neuron. The results suggest that the set of conductances expressed in the various neuron types may be divided into two groups: one group is responsible for the major characteristics of the firing behavior in each population and the other is more involved with degeneracy. These models provide experimentally testable predictions on the combination and relative proportion of the different conductance types that should be present in hippocampal CA1 pyramidal cells and interneurons.
keywords: medicine and health sciences
  action potentials
  membrane potential
  electrophysiology
  neuroscience
  optimization
  mathematics
  ganglion cells
  interneurons
  neuronal dendrites
  animal cells
  cellular neuroscience
  neuronal morphology
  cell biology
  pyramidal cells
  physiology
  neurons
  biology and life sciences
  cellular types
  physical sciences
  neurophysiology
title: The physiological variability of channel density in hippocampal CA1 pyramidal cells and interneurons explored using a unified data-driven modeling workflow