A multi-compartment model of a thalamic interneuron implemented in NEURON¶
In this example we illustrate how Uncertainpy can be used on models implemented
in NEURON.
For this example, we select a previously published model of an interneuron in
the dorsal lateral geniculate nucleus Halnes et al., 2011.
Since the model is in implemented in NEURON,
the original model can be used directly with Uncertainpy with the use of
NeuronModel.
The code for this case study is found in
/examples/interneuron/uq_interneuron.py
.
To be able to run this example you require both the NEURON simulator,
as well as the interneuron model saved in the folder /interneuron_model/
.
In the original modeling study, a set of 7 parameters were tuned manually through a series of trials and errors until the interneuron model obtained the desired response characteristics. The final parameter set is:
Parameter | Value | Unit | Neuron variable | Meaning |
---|---|---|---|---|
\(g_{\mathrm{Na}}\) | 0.09 | \(\text{S/cm}^2\) | gna |
Max \(\text{Na}^+\)-conductance in soma |
\(g_{\mathrm{Kdr}}\) | 0.37 | \(\text{S/cm}^2\) | gkdr |
Max direct rectifying \(\text{K}^+\)-conductance in soma |
\(g_{\mathrm{CaT}}\) | 1.17e-5 | \(\text{S/cm}^2\) | gcat |
Max T-type \(\text{Ca}^{2+}\)-conductance in soma |
\(g_{\mathrm{CaL}}\) | 9e-4 | \(\text{S/cm}^2\) | gcal |
Max L-type \(\text{Ca}^{2+}\)-conductance in soma |
\(g_{\mathrm{h}}\) | 1.1e-4 | \(\text{S/cm}^2\) | ghbar |
Max conductance of a non-specific hyperpolarization activated cation channel in soma |
\(g_{\mathrm{AHP}}\) | 6.4e-5 | \(\text{S/cm}^2\) | gahp |
Max afterhyperpolarizing \(\text{K}^+\)-conductance in soma |
\(g_{\mathrm{CAN}}\) | 2e-8 | \(\text{S/cm}^2\) | gcanbar |
Max conductance of a \(\text{Ca}^{2+}\)-activated non-specific cation channel in soma |
To perform an uncertainty quantification and sensitivity analysis of this model, we assume each of these 7 parameters have a uniform uncertainty distribution in the interval \(\pm 10\%\) around their original value. We create these parameters similar to how we did in the Hodgkin-Huxley example:
# Define a parameter list
parameters= {"gna": 0.09,
"gkdr": 0.37,
"gcat": 1.17e-5,
"gcal": 0.0009,
"ghbar": 0.00011,
"gahp": 6.4e-5,
"gcanbar": 2e-8}
# Create the parameters
parameters = un.Parameters(parameters)
# Set all parameters to have a uniform distribution
# within a 20% interval around their fixed value
parameters.set_all_distributions(un.uniform(0.2))
A point-to-point comparison of voltage traces is often uninformative, and we therefore want to perform a feature based analysis of the model. Since we examine a spiking neuron model, we choose the features in SpikingFeatures:
# Initialize the features
features = un.SpikingFeatures(features_to_run="all")
We study the response of the interneuron to a somatic current injection
between \(1000 \text{ ms} < t < 1900 \text{ ms}\).
SpikingFeatures
needs to know the start and end time of this
stimulus to be able to calculate certain features.
They are specified through the stimulus_start
and
stimulus_end
arguments when initializing NeuronModel
.
Additionally, the interneuron model uses adaptive time steps,
meaning we have to set interpolate=True
.
In this way we tell Uncertainpy to perform an interpolation to get the
output on a regular form before performing the analysis:
We also give the path to the folder where
the neuron model is stored with path="interneuron_model/"
.
NeuronModel
loads the NEURON model from mosinit.hoc
,
sets the parameters of the model,
evaluates the model and returns the somatic membrane potential of the neuron,
(the voltage of the section named "soma"
).
NeuronModel
therefore does not require a model function.
# Initialize the model with the start and end time of the stimulus
model = un.NeuronModel(path="interneuron_model/", interpolate=True,
stimulus_start=1000, stimulus_end=1900)
We set up the problem, adding our features before we use polynomial chaos expansion with point collocation to compute the statistical metrics for the model output and all features. We also set the seed to easier be able to reproduce the result.
# Perform the uncertainty quantification
UQ = un.UncertaintyQuantification(model,
parameters=parameters,
features=features)
# We set the seed to easier be able to reproduce the result
data = UQ.quantify(seed=10)
The complete code becomes:
import uncertainpy as un
# Define a parameter list
parameters= {"gna": 0.09,
"gkdr": 0.37,
"gcat": 1.17e-5,
"gcal": 0.0009,
"ghbar": 0.00011,
"gahp": 6.4e-5,
"gcanbar": 2e-8}
# Create the parameters
parameters = un.Parameters(parameters)
# Set all parameters to have a uniform distribution
# within a 20% interval around their fixed value
parameters.set_all_distributions(un.uniform(0.2))
# Initialize the features
features = un.SpikingFeatures(features_to_run="all")
# Initialize the model with the start and end time of the stimulus
model = un.NeuronModel(path="interneuron_model/", interpolate=True,
stimulus_start=1000, stimulus_end=1900)
# Perform the uncertainty quantification
UQ = un.UncertaintyQuantification(model,
parameters=parameters,
features=features)
# We set the seed to easier be able to reproduce the result
data = UQ.quantify(seed=10)