My own preference would be to make Erev be the dependent variable
and that is what I have done in this implementation of an event-driven synaptic mechanism based on the NMODL code in your first post. I have inserted explanatory comments where this implementation differs from your original code.
Code: Select all
NEURON {
POINT_PROCESS gabaa
USEION cl READ ecl WRITE icl VALENCE -1
NONSPECIFIC_CURRENT ihco3
: value of C is now from element 0 of the NetCon's weight vector
: POINTER C
: replace the idiosyncratic curr with the widely used i
: STATEs are automatically range, so no need to declare this
: RANGE C0, C1, C2, D1, D2, O1,O2, curr, hco3i, hco3o
RANGE i, hco3i, hco3o
: the WRITE icl declaration in the USEION statement
: does not make it possible to discover the value of the chloride
: current generated by any particular instance of this mechanism.
: make this value accessible to hoc and Python by declaring icl RANGE
RANGE icl
: synaptic delay parameter del is now handled
: by the NetCon that delivers events to this mechanism
: RANGE g, gmax, rb1, rb2, q0, tr, del, tdur
RANGE g, gmax, rb1, rb2, q0, tr, tdur
: the values of ecl, ehco3, Pr (the fraction of synaptic conductance
: that is attributable to hco3), and Erev are not independent of each
: other, because
: Erev = Pr * ehco3 + (1 - Pr) * ecl
: furthermore, since the value of Pr is constrained to lie in the
: range 0..1, it follows that Erev is limited to the range ehco3..ecl
: Erev is most easily measured in an experiment, but its value
: is a consequence of other values that have direct biological bases
: so it becomes an ASSIGNED
: also, its value might easily be a function of spatial location,
: so it should be declared RANGE
: GLOBAL Erev
: finally, the value of ehco3 could easily be of interest
: so make it RANGE too
RANGE Erev, ehco3
GLOBAL Rb1, Rb2, Ru1, Ru2, Rd1, Rd2, Rr1, Rr2, a1, b1, a2, b2, p, q
}
UNITS {
(nA) = (nanoamp)
(mV) = (millivolt)
(pS) = (picosiemens)
(umho) = (micromho)
(mM) = (milli/liter)
(uM) = (micro/liter)
(mA) = (milliamp)
F = (faraday) (coulomb)
R = (mole k) (mV-coulomb/degC)
}
PARAMETER {
: VERIFY THE FOLLOWING!
: NEURON will ignore this assignment because in NEURON
: celsius is a top-level global parameter that is accessible
: to hoc and Python as celsius and h.celsius respectively
: celsius = 37 (degC) : this value will be ignored
: Erev is now an ASSIGNED with value calculated from ecl, ehco3, and Pr
: see comment in NEURON block
: Erev = -80 (mV) : reversal potential
hco3i = 16 (mM)
hco3o = 26 (mM)
gmax = 28 (pS) : maximal conductance
tdur = 5 (ms)
: synaptic delay is specified by the NetCon that delivers the event
: that activates this synaptic mechanism
: del = 0 (ms)
Pr = 0.3
dmax = 3.6 (/ms)
kh = 2 (mM)
nh = 2
: Rates
rb1 = 30 (/mM /ms)
Ru1 = 90e-03 (/ms)
rb2 = 15 (/mM /ms)
Ru2 = 180e-03 (/ms)
b1 = 0.2 (/ms)
a1 = 1.1 (/ms)
b2 = 2.5 (/ms)
a2 = 0.142 (/ms)
Rd1 = 0.006 (/ms)
Rr1 = 0.35e-03 (/ms)
Rr2 = 0.05 (/ms)
p = 0.001 (/ms)
q0 = 1e-08 (/mM /ms)
}
ASSIGNED {
v (mV) : postsynaptic voltage
: use i instead of the idisyncratic curr
: curr (nA) : current = g*(v - Erev)
i (nA) : total current generated by an instance of this mechanism
g (pS) : conductance
tr (mM)
ecl (mV)
ehco3 (mV)
Rb1 (/ms)
Rb2 (/ms)
Rd2 (/ms)
q (/ms)
icl (nA)
C (mM)
ihco3 (nA)
: a state variable that reflects the presence or absence
: of transmitter in the synaptic gap
on (1) : 1 if synapse has been activated, 0 otherwise
celsius (degC) : see comment in PARAMETER block
: see comment in NEURON block
Erev (mV) : apparent synaptic reversal potential
: calculated from ecl, ehco3, and Pr
}
STATE {
C0 : unbound
C1 : single GABA bound
C2 : double GABA bound
O1 : Open state from single bound state
O2 : Open state from single bound state
D1 : single GABA bound, desensitized
D2 : double GABA bound, desensitized
}
INITIAL {
C0 = 1
C1 = 0
C2 = 0
O1 = 0
O2 = 0
D1 = 0
D2 = 0
tr = 0
on = 0
: if this mechanism is used in a model that assumes hco3i and hco3o
: are constant during a simulation, just calculate ehco3 here
ehco3 = log(hco3i/hco3o)*(celsius + 273.15)*R/F
Erev = Pr * ehco3 + (1 - Pr) * ecl : see comment in NEURON block
}
BREAKPOINT {
SOLVE kstates METHOD sparse
: see comment in INITIAL block
: ehco3 = log(hco3i/hco3o)*(celsius + 273.15)*R/F
g = gmax * (O1 + O2)
: use i instead; see comment in NEURON block
: curr = (1e-6) * g * (v - Erev)
: these are incorrect unless ehco3 == Erev by luck or contrivance
: ihco3 = Pr*curr
: icl = (1-Pr)*curr
: to compute i, calculate the individual current components
: and add them together
: these are always correct
icl = (1e-6) * g * (1-Pr) * (v - ecl)
ihco3 = (1e-6) * g * Pr * (v - ehco3)
i = icl + ihco3
}
KINETIC kstates {
: FUNCTION release() has been replaced by statements in the
: NET_RECEIVE block--see NET_RECEIVE below
: release()
Rb1 = rb1 * tr
Rb2 = rb2 * tr
Rd2 = dmax/( 1 + (kh/tr)^nh)
q = q0 * tr
~ C0 <-> C1 (Rb1, Ru1)
~ C1 <-> C2 (Rb2, Ru2)
~ C1 <-> O1 (b1, a1)
~ C2 <-> O2 (b2, a2)
~ O1 <-> D1 (Rd1, Rr1)
~ O2 <-> D2 (Rd2, Rr2)
~ D1 <-> D2 (q, p)
CONSERVE C0+C1+C2+D1+D2+O1+O2 = 1
}
COMMENT
replaced by statements in the NET RECEIVE block
FUNCTION release(){
at_time(del)
at_time(del+tdur)
if (t>del && t<=(del+tdur)) {
tr = C
}
else{ tr = 0}
}
ENDCOMMENT
NET_RECEIVE (C (mM), on) {
if (flag == 0) { : this event was triggered by
: a presynaptic variable crossing a NetCon's threshold
: in a positive-going direction (e.g. by a presynaptic spike)
if (!on) { : "activate" the synapse
tr = C : concentration of transmitter in synaptic gap equals C
net_send(tdur, 1) : a self-event that will signal when
: tr should return to zero
on = 1 : so it is clear that the synapse has been activated
} else {
: if the synapse is already active,
: move the self-event so that the transmitter persists
: in the synaptic gap until tdur ms after the most recent
: event with flag==0
net_move(t + tdur)
}
}
if (flag == 1) {
: the self-event has arrived
: time for tr to fall back to 0
tr = 0
on = 0
}
: given the complexity of the gating mechanism
: superposition can't be exploited--
: a separate instance of this mechanism
: will be needed for each afferent event stream
}
Some comments about the performance of this mechanism with its present parameter values:
1. ehco3 is about -11.7 mV. ecl, however, will be very close to the model cell's resting potential. Consequently ihco3 will be much larger than icl, and synaptic activation will depolarize the cell.
That said, a single synaptic activation has little effect on membrane potential unless transmitter concentration in the synaptic gap is fairly large. Given a single compartment model cell that has surface area 100um2 and the pas mechanism with g_pas 0.001 S/cm2 and e_pas -70 mV (no other channels, so resting potential is -70 mV), a single synaptic activation that raises transmitter concentration in the synaptic gap to 0.1 mM (NetCon's weight[0] equal to 0.1) for 5 ms depolarizes the cell by only about 0.43 mV.
2. A single synaptic activation that raises transmitter concentration tr in the synaptic gap to 1 mM for 5 ms elicits a synaptic conductance transient with a biphasic time course consisting of a brief early peak followed by a later, very broad peak of similar amplitude. The amplitudes and times of these peaks are 16 pS at 1.85 ms after the onset of the increase of tr, and 18 pS at about 33.175 ms after the increase of tr (about 27 ms after transmitter vanishes from the gap). They are separated by a deep minimum of about 2.85 pS at about 6 ms after the increase of tr (about 1 ms after transmitter vanishes from the gap). Is there any experimental evidence for such a time course at any GABA-A synapse?
I'll be glad to provide by email the source code for the model that produced these results, if you like.
My guess is that the bicarbonate concentrations aren't quite right.