time shift in SEClamp

The basics of how to develop, test, and use models.
Post Reply
cafischer

time shift in SEClamp

Post by cafischer »

I played a vector into a SEClamp and then recorded the membrane potential to see whether it matches the played vector but it did not. Firstly the beginning can be quite different if v_init does not match the first value in the played vector. But also the second value still did not match. Secondly, all the following values resembled the played vector more but still deviated by a certain amount.

Why is this the case?
How can one get exactly the membrane potential that was feeded into the SEClamp?
In the end I also want the current traces to correspond to the current produced at the clamped potential. How can the current trace be matched to the clamped potential? (probably the same answer as for the question above).

The following code demonstrates the problem:

Code: Select all

from neuron import h
import numpy as np
import pylab as pl
h.load_file("stdrun.hoc")  # load NEURON libraries
h("""cvode.active(0)""")  # unvariable time step in NEURON

# create simple cell model
cell = h.Section(name='soma')
cell.cm = 1
cell.Ra = 100
cell.diam = 100
cell.L = 100
cell.nseg = 1
cell.insert('pas')

# create membrane potential trace
tstop = 50
dt = 0.025
t = np.arange(0, tstop+dt, dt)
v = np.linspace(-80, -50, len(t))

t_clamp = h.Vector()
t_clamp.from_python(t)
v_clamp = h.Vector()
v_clamp.from_python(v)

# seclamp
stim = h.SEClamp(0.5, sec=cell)
stim.rs = 1e-3  # series resistance should be much smaller than input resistance of the cell
stim.dur1 = 1e9
v_clamp.play(stim._ref_amp1, dt)

# recordings
ipas = h.Vector()
ipas.record(cell(.5).pas._ref_i)

v_rec = h.Vector()
v_rec.record(cell(.5)._ref_v)

t_rec = h.Vector()
t_rec.record(h._ref_t)

# run simulation
h.tstop = t[-1] + dt
h.dt = dt
# h.v_init = np.array(v_clamp)[0]
h.init()
h.run()

# plot
f, (ax1, ax2) = pl.subplots(2, 1, sharex=True)
ax1.plot(np.array(t_rec), np.array(v_rec), 'k', label='recorded v')
ax1.plot(np.array(t_clamp), np.array(v_clamp), 'b', label='clamped v')
ax1.set_ylabel('Membrane\npotential (mV)')
ax2.plot(np.array(t_rec), np.array(ipas), 'k', label='ipas')
ax2.set_ylabel('Current (mA/cm2)')
ax2.set_xlabel('Time (ms)')
ax1.legend()
pl.show()
ted
Site Admin
Posts: 6286
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: time shift in SEClamp

Post by ted »

Hint: initialization and series resistance. If you're a physicist or electrical engineer, those two clues should be enough to get you on the right track.

If not--
cafischer wrote:the beginning can be quite different if v_init does not match the first value in the played vector.
v_init specifies membrane potential at t=0, and (indirectly) the amount of charge stored in membrane capacitance. Rate of change of membrane potential is proportional to capacitive current. If the membrane has no ion channels, capacitive current equals clamp current, and clamp current equals
(clamp's command potential - membrane potential)/clamp's series resistance
Numerator is always finite, denominator is never 0, so membrane potential always lags behind clamp's command potential.
How can one get exactly the membrane potential that was feeded into the SEClamp?
Block all ion channels, and reduce both the clamp's series resistance and the membrane capacitance to 0. Not going to happen on this planet. You can, however, make the series resistance very small--read the documentation on the SEClamp class to discover the name of this parameter and its default value. And you can also reduce all ionic conductances to 0. You could reduce specific membrane capacitance to a very small value, but then you'd be dealing with very nonphysiological values.
In the end I also want the current traces to correspond to the current produced at the clamped potential.
By "current traces", which current(s) do you mean? And what is "the current produced at the clamped potential"?
cafischer

Re: time shift in SEClamp

Post by cafischer »

As I understand it know without ion channels the equation is
cm * dV/dt = i_clamp = (V_clamp - Vm) / R_clamp

Assuming you use the forward-Euler (V_new = V_old + dV/ dt * dt) for the integration (which is not the case in NEURON) then:
V_new = V_old + (V_clamp - V_old) * dt / (cm * R_clamp)
= dt/(cm*R_clamp) * V_clamp + (1 - dt/(cm*R_clamp)) * V_old
which means that V_new = V_clamp(old) if dt/cm*R_clamp = 1. So that the Vm trace is equal to the V_clamp trace if we shift it by one time step.

If you add ionic currents the equation becomes more complicated:
cm * dV/dt = -sum i_ion + (V_clamp - Vm) / R_clamp
so that
V_new = -1/cm * sum i_ion + dt/(cm*R_clamp) * V_clamp + (1 - dt/(cm*R_clamp)) * V_old

So the voltage will be driven away from the clamped potential by the ionic currents. But if you want to eliminate the influence of V_old wouldn't I still want dt/(cm*R_clamp) = 1 implying R_clamp = dt/cm?
By the way I did not find a default value of the series resistance in the documentation but printing it gave me 1MOhm.

Why I look at this is because I want to know the current flowing through a certain ion channel inserted in the membrane in response to a given membrane potential ideally resembling the clamped potential.
So as I know see it I get the closest estimate if: R_clamp = dt/cm and if I shift the membrane potential and current trace by 1 time step.
Would you say that this is correct?
ted
Site Admin
Posts: 6286
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: time shift in SEClamp

Post by ted »

I did not find a default value of the series resistance in the documentation but printing it gave me 1MOhm.
Good.
As I understand it . . .
. . . wouldn't I still want dt/(cm*R_clamp) = 1 implying R_clamp = dt/cm?
. . . So as I know see it I get the closest estimate if: R_clamp = dt/cm and if I shift the membrane potential and current trace by 1 time step.
No. You're making things more complicated than necessary, and reaching incorrect conclusions, including completely unsupported inferences about what integration time step to choose and mystical ideas about why it is or isn't necessary to shift this or that recording of some current or other by one time step. Throw all that out and start with a blank page.

First, the physics of the situation dictate that in the actual physical system there will always be a lag between the time course of injected current and the resulting change of membrane potential. This is because membrane potential can change only if charge is removed from, or deposited on, membrane capacitance. This must take a finite amount of time unless you are able to apply an infinitely large current (which is impossible).

Second, finite difference numerical integration inherently results in a lag between the time course of injected current and the resulting change of membrane potential in the computational model. To see why, write the differential equation with the derivative on the left, and everything else on the right, so that it is of the form
dv/dt = f(v,t)
To do that, rearrange this
cm*dv/dt + i_ionic = injected_current
so that it becomes
dv/dt = (injected_current - i_ionic)/cm

For a model with linear passive ion channels, this is
dv/dt = (injected_current - g_pas*(v - e_pas))/cm

Let's keep things simple by assuming that the injected current is generated by a current clamp. Then the injected current is independent of v and depends only on time t, and the equation becomes
dv/dt = (i(t) - g_pas*(v(t) - e_pas))/cm
We can simplify this further by assuming that e_pas is 0
Equation 1: dv/dt = (i(t) - g_pas*v(t))/cm

The explicit Euler method approximates the solution to dv/dt = f(v,t) by
v(t+Dt) = v(t) + Dt*f(v(t),t)
where Dt is just a lazy way of writing "delta t".

Applying explicit Euler to Equation 1 yields this recurrence equation
Equation 2: v(t+Dt) = v(t) + (i(t) - g_pas*v(t))*Dt/cm

By inspection it is clear that the injected current at the present moment does not have an immediate effect on membrane potential--i(t) can only affect v at future times.

The implicit, or "backward", Euler method approximates the solution to dv/dt = f(v,t) by
v(t+Dt) = v(t) + f(v(t+Dt),t+Dt)*Dt
Applying implicit Euler to Equation 1, after some rearrangement and grouping of terms, gives us this recurrence equation
Equation 3: v(t+Dt) = v(t)/(1 + Dt*g_pas/cm) + i(t+Dt)*Dt/cm/(1 + Dt*g_pas/cm)

Once again it is clear that the injected current at the present moment does not have an immediate effect on membrane potential--i(t) can only affect v at future times.

So there is always a lag between the injected current and the effect on membrane potential. Furthermore, this is the case even if e_pas is nonzero, and it is still the case if the driving function is current delivered by a voltage clamp that has finite series resistance.
I want to know the current flowing through a certain ion channel inserted in the membrane in response to a given membrane potential ideally resembling the clamped potential.
Easy to do in a computational model. Vector play allows one to force an SEClamp's command potential to follow a piecewise linear approximation to an arbitrary time course. Just make sure that the SEClamp's series resistance is much smaller than the input resistance of the model cell, so that local voltage control is good. As long as the NMODL code that specifies the voltage-gated ion channel's properties exposes the current through those channels to hoc or Python, you're good to go. If it doesn't, that can be easily fixed. If the channel's properties are specified by an instance of the Channel Builder class, nothing special will need to be done.
cafischer

Re: time shift in SEClamp

Post by cafischer »

Thanks for the elaborate explanations!
I was thinking about the last part again. I said that:
I want to know the current flowing through a certain ion channel inserted in the membrane in response to a given membrane potential ideally resembling the clamped potential.
As you explained that the SEClamp always will have a time lag I wondered whether it is possible to feed the voltage directly into the channel model and compute the error so that this source of error is eliminated. Is it possible to do this in NEURON? If yes, how would you do it?
ted
Site Admin
Posts: 6286
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: time shift in SEClamp

Post by ted »

On the topic of series resistance--you might find it useful to get and read The Axon Guide from Molecular Devices (free download). Also they have a webinar with the promising title "Series Resistance Compensated or Not" (haven't seen it myself, so can't vouch for its contents).
I want to know the current flowing through a certain ion channel inserted in the membrane in response to a given membrane potential ideally resembling the clamped potential.
As you explained that the SEClamp always will have a time lag I wondered whether it is possible to feed the voltage directly into the channel model and compute the error so that this source of error is eliminated. Is it possible to do this in NEURON? If yes, how would you do it?
There's always going to be a difference between the command potential and the local membrane potential where the patch electrode is attached. Just make the SEClamp's series resistance very small, and the difference will be negligible.
Post Reply