Explanation of the Matlab functions in the stocHHastic package
The attached Matlab code implements the stochastic HodgkinHuxley
model with ionchannel gating modeled as Markov chains. We provide
both the full Markov chain model as well as its stochasticshielding
approximation (folder HH). In addition, we provide the code comparing
the theoretical covariance matrices for the transitions and occupation
numbers in the full and simplified models with the covariance matrices
obtained from the simulations (folder QandC).
1.Folder HH
Programs with "full_mark" in the title utilize a full markov process
for ion channel stochasticity; "simp_mark" in the title have the
shielding approximation.
Files with "_Vclamp" at the end are entirely in voltage clamp and must
have the voltage given to the program, while files with the "_Iclamp"
suffix are in current clamp. Please see the comments in the
files for details of parameter values.
The eight common inputs parameters to both current clamp and voltage
clamp functions are:
dt: timestep in ms,
T: total time interval (duration of simulated signal) in ms,
Nanoise and Knoise: binary variables representing whether or not
sodium and potassium, respectively, are stochastic or deterministic,
Nanum and Knum: channel densities for sodium and potassium per m2,
trans: duration of the transient (initial deterministic period that IS
still included in the output),
area: membrane patch size in m2.
Scripts with "_Vclamp" at the end also have a ninth input,
v: the voltage at which they are clamped in mV.
The current clamp scripts return a value for the membrane potential
that the voltage clamp scripts do not, but after that both output the
same nine values; total current, sodium current, potassium current, a
timetrack vector that is in seconds, a sodium matrix showing the
number of channels in each Markov state, a potassium matrix showing
the number of channels in each Markov state, the total number of
sodium channels, the total number of potassium channels, and the time
it took the simulation to run.
Figure 2 is in voltage clamp, and was created with the respective
files with the suffix "_Vclamp" using the parameters given in the
figure caption. Figure 3 was created with the Matlab files with suffix
"_Iclamp" using parameters in the figure caption.
2.Expanded Model
These files are very similar in notation and style to the HH models,
with "full_mark" and "simp_mark" representing programs without and
with the shielding approximation, respectively.
All of the input and output parameters described for the HH package
are the same, however there are additional inputs
Pnum and SKnum: channel densities for P/Q channels and SK
channels, respectively
Idc: Input current (should be a vector of equal length as
there are timesteps)
g_input: Input conductance (with reversal potential of 0mV)
Buffering_constant: The buffering constant represents a
constant that the amount of calcium entering the cell is
divided by, to account for the endogenous buffering of the
cell. Since the calcium concentration will only be used to
gate the SK channels, this is a valid approximation.
Catau: The time constant at which calcium is removed from the
cell, modeled as an exponential decay.
And additional outputs
I_pca and I_sk: P/Q calcium and SK currents
Caconc: the concentration of calcium in the cell
PCa and SK: the number of channels in each Markov State
throughout the trajectory
g_current: the current due to the g_input conductance
Spiketriggered averages are computed as the normalized
crosscorrelogram of the input current with the spike events
(represented as delta functions).
To recreate responses to current or conductance input, simply run the
files "OUnoise" for OrnsteinUhlenbeck noise or "poisson_input" for
Poisson conductance events (note that OUnoise needs to be called as a
function but poisson_input is simply run as a terminal script), and
use the outputs as the Idc or g_input input for the model script. The
other variable should be simply set to a vector of zeros.
Note that the HodgkinHuxley Files in the other package do not have an
input current value, so to recreate those results, simply set the
values for Pnum and SKnum to 0 in the expanded model
3.Folder QandC
The mfile names are selfexplanatory. Files with prefix "simulations"
compute the covariance matrices Q and C from the simulated data of
stochastic ion channel gating in voltage clamp, whereas files with
prefix "theoretical" compute those matrices from the theoretical
expressions shown in the manuscript. File names containing "full_MC"
use the full Markov chain model and files names containing "simp_MC"
use the stochasticshielding approximation. Although three of the four
mfiles are defined as functions, they do not take input parameters,
as they are defined in the within the file. These files where defined
as functions so that they can call functions defined within the files.
The program "theory_vs_sims_CovN_MC.m" compares the results from
theory and simulations using the output files (matfiles) of the other
three mfiles. The matfiles used for the manuscript are provided as
an example, along with .eps figures showing the results.
By Nicolaus T. Schmandt (n.schmandt@gmail.com) & Roberto Fernandez
Galan (rfgalan@case.edu),
Case Western Reserve University, July 2012.
