www.neuron.yale.edu

I am using the phase plane plot tool in NEURON GUI. I find if I want to plot dV/dt ~ V, I can plot total current of that segment versus V. However, is there a variable which I can access from hoc that is exactly this dV/dt value? If that is possible, it would be much easier for me to plot also, say, dg/dt value.

Otherwise I would have to record the vectors down and make difference of the vectors, which might be troublesome especially when I am using variable step integration in my simulation.

eacheon wrote:I am using the phase plane plot tool in NEURON GUI. I find if I want to plot dV/dt ~ V, I can plot total current of that segment versus V.

You are of course referring to total ionic current, not total membrane current.

However, is there a variable which I can access from hoc that is exactly this dV/dt value?

Numerical integration always produces approximate results. The best one can get is results
that are internally consistent and within an acceptable error tolerance.

NEURON's default integrator is implicit Euler, which has first order accuracy. Greater
precision is achievable by setting secondorder = 2, which gives second order accurate
results for voltage and currents--but at different times in the course of the simulation!
That is, v is 2nd order correct at times n*dt, and i is 2nd order correct at n*dt - dt/2
(see
http://www.neuron.yale.edu/neuron/stati ... econdorder
). If you must have both correct to the highest precision at the same time, use adaptive
integration (aka cvode).

Otherwise I would have to record the vectors down and make difference of the vectors, which might be troublesome especially when I am using variable step integration in my simulation.

Not at all. Use
vdest.record(&var, DT)
where DT is your desired sampling interval, and NEURON will generate the solution wth
adaptive integration but forcing the integrator to stop at integer multiples of DT to calculate
and record the variable to your vector. So you have evenly spaced samples, and can
readily use the Vector class's deriv() method to calculate your derivative.

The only caveat is to use cvode.record rather than vector.record() if you are using local
variable time step integration.

ted wrote: You are of course referring to total ionic current, not total membrane current.

Exactly.

ted wrote: Numerical integration always produces approximate results. The best one can get is results
that are internally consistent and within an acceptable error tolerance.

NEURON's default integrator is implicit Euler, which has first order accuracy. Greater
precision is achievable by setting secondorder = 2, which gives second order accurate
results for voltage and currents--but at different times in the course of the simulation!
That is, v is 2nd order correct at times n*dt, and i is 2nd order correct at n*dt - dt/2
(see
http://www.neuron.yale.edu/neuron/stati ... econdorder
). If you must have both correct to the highest precision at the same time, use adaptive
integration (aka cvode).

Thanks for the very patient and detailed answer on this. I've read similar notes provided by you on this forum in two or three other threads.

So you are saying there's no way I can access the dg/dt values directly from hoc? It is a pity since anyway NEURON's integrator knows these values at a first hand.

ted wrote: Not at all. Use
vdest.record(&var, DT)
where DT is your desired sampling interval, and NEURON will generate the solution wth
adaptive integration but forcing the integrator to stop at integer multiples of DT to calculate
and record the variable to your vector. So you have evenly spaced samples, and can
readily use the Vector class's deriv() method to calculate your derivative.

The only caveat is to use cvode.record rather than vector.record() if you are using local
variable time step integration.

Thanks for this nice tip! I will do as you suggested. This way I may actually directly using the "phase plane" tools to plot those in NEURON...

So you are saying there's no way I can access the dg/dt values directly from hoc?

Correct. Whether dg/dt is even computed internally depends on how
the channel's properties have been specified, and on which integrator
is being used.

I will do as you suggested. This way I may actually directly using the "phase plane" tools to plot those in NEURON...

You will have to plot your own phase plane using the Vector class's
plot() method, with the syntax vec.plot(graphobj, x_vec) (see
https://www.neuron.yale.edu/neuron/stat ... .html#plot
), not with the phase plane graph created by
NEURON Main Menu / Graph / Phase Plane

eacheon wrote:

ted wrote: You are of course referring to total ionic current, not total membrane current.

Exactly.

Ted, I might miss something here. Shouldn't I be plotting i_cap versus v as the dV/dt ~ V phase plot? And I am thinking i_cap (for a segment in a section) almost always does not equal to total ionic current (for that segment), am I right, or I am missing some serious thing here?

eacheon wrote:Ted, I might miss something here. Shouldn't I be plotting i_cap versus v as the dV/dt ~ V phase plot?

It looks like we were both missing something. dV/dt = membrane capacitive
current divided by membrane capacitance. Total membrane current is the
sum of ionic and capacitive components, so in a model that has only one
compartment, dV/dt = - ionic membrane current divided by membrane
capacitance. But in a model with more than one compartment, there are
axial currents that destroy this nice relationship, i.e. as you put it

i_cap (for a segment in a section) almost always does not equal to total ionic current (for that segment)

So the strategy that works with all models is to compute dV/dt numerically.

The deep question is: what sense can one make out of a phase plane plot
in a model that has multiple compartments?

I seem to have sowed confusion. These statements are true:

dV/dt = membrane capacitive
current divided by membrane capacitance. Total membrane current is the
sum of ionic and capacitive components, so in a model that has only one
compartment, dV/dt = - ionic membrane current divided by membrane
capacitance. But in a model with more than one compartment, there are
axial currents that destroy this nice relationship, i.e. as you put it

i_cap (for a segment in a section) almost always does not equal to total ionic current (for that segment)

But the range variable i_cap _is_ the current through the membrane capacitance, so the ratio
i_cap(x)/(area(x)*cm(x)) at any internal node
should be equal to dv/dt at that node. Numerical differentiation of v is not necessary.

I am still not 100% sure what icap is. I made a test with a single-compartment model where I compared icap with the negative sum of all ionic currents plus the injected current and found that there where differences between those two which reduced with increasing step size. This suggests that icap is not just

the current through the membrane capacitance

but is computed by numerical differentiation of the integrated membrane potential V times the cells capacitance.

Is this true?
If not, why is it not equal to the sum of the currents?
Is there a variable in NEURON from which I can directly access dV/dt before integration of V?

cafischer wrote:I am still not 100% sure what icap is. I made a test with a single-compartment model where I compared icap with the negative sum of all ionic currents plus the injected current and found that there where differences between those two which reduced with increasing step size. This suggests that icap is not just

the current through the membrane capacitance

but is computed by numerical differentiation of the integrated membrane potential V times the cells capacitance.

Is this true?

Not at all. If the model is passive and is subjected to a tonic (static, unchanging) stimulus, that is exactly what should happen with implicit Euler. It has no implications for how i_cap is computed.

The fixed step method implementation does not expose the derivatives of the states. But they are available with the variable step integrators. You need to consider

http://neuron.yale.edu/neuron/static/ne ... de.dstates

also helpful in this context is

http://neuron.yale.edu/neuron/static/ne ... .statename