Time dependent conductance in DERIVATIVE block

NMODL and the Channel Builder.
Post Reply
myreil
Posts: 11
Joined: Mon Mar 14, 2011 2:25 pm

Time dependent conductance in DERIVATIVE block

Post by myreil »

Hi!
I am new in NEURON and want to add a channel in a hh model.
My problem is that in the equations in the DERIVATIVE block the variable change depending on time. how can I define different values to those variables as time course passes?
Something like:

Code: Select all

if (t < thres) {
		Ga2=r2*(1-exp(-t/tau))
	} else{
		Ga2=r2*((exp(-t-thres)/tau)-exp(-t/tau))
	}
}
Values of Ga1 and Ga2 are used to solve the differential equations in the DERIVATIVE block

Code: Select all

DERIVATIVE states {
	O1'=Ga1-(Gd1+e_ct+Ga1)*O1+(e_ct-Ga1)*O2-Ga1*C1
	O2'=e_ct*O1-(e_tc+Gd2)*O2+Ga2*C1
	C1'=Gd2*O2-(Gr+Ga2)*C1
}
Hope that this makes sense to you.
Thank you!
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Time dependent conductance in DERIVATIVE block

Post by ted »

Not much, but it's Monday.

What is the conceptual model of this mechanism? If it was originally stated as a family of states with transitions between them, it might be easier to use a KINETIC block rather than to try to rewrite it in the form of a set of ODEs.
myreil
Posts: 11
Joined: Mon Mar 14, 2011 2:25 pm

Re: Time dependent conductance in DERIVATIVE block

Post by myreil »

Dear Ted,
Thank you for your answer.

I just want to change the values of Ga1 and Ga2 as time passes.
Like:
Ga1= something, when ton<t<toff
Ga1=something else, when t>toff

My question is if I can have time as a threshold to change values? Because most examples I have seen use voltage or temperature.

And I have only the ODEs, I think a kinetic model would be more complicated.
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Time dependent conductance in DERIVATIVE block

Post by ted »

Since you have the analytical description of Ga2, and the value of Ga2 is used in just one equation, simply add a FUNCTION block to your mod file

Code: Select all

FUNCTION Ga2() {
  if ( . . . ) {
    . . .
  }
}
You will probably also have to declare
t (ms)
in the ASSIGNED block. The resulting code will run fine with fixed time step integration, but will not work reliably with adaptive integration.

If you want something that works reliably with adaptive integration, let me know and I'll show you how.
myreil
Posts: 11
Joined: Mon Mar 14, 2011 2:25 pm

Re: Time dependent conductance in DERIVATIVE block

Post by myreil »

Thank you!! It works like that! If I use FUNCTION t needs to be declared in the ASSIGNED block, but I did it also with PROCEDURE and in this case t cannot be declared in the ASSIGNED block.

If you want something that works reliably with adaptive integration, let me know and I'll show you how.[/quote]

Yes,please. That would be of great help.
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Time dependent conductance in DERIVATIVE block

Post by ted »

In brief outline, here is the problem, and a clue to what must be done to make your code work with adaptive integration.

Adaptive integration may change the time step in the middle of a simulation. If you have something that must happen at a particular time tthresh, you must ensure that it really does happen at that time--because the adaptive integrator might increment t by a value that steps well past tthresh.
if (t>=tthresh)
will almost certainly fail to detect the exact instant at which t==tthresh.

The solution is to use events. In the INITIAL block insert a net_send(tthresh, 1) statement. This produces a self-event that returns with flag==1 when t==tthresh. The mechanism will also need a NET_RECEIVE block that tells it what to do when this self-event arrives. For example, if this mechanism has a PARAMETER called a, and we want a to become 1 when the self-event arrives, this NET_RECEIVE block would do the trick:

Code: Select all

NET_RECEIVE (w) {
  if (flag==1) { : is this the self-event that we sent at initialization?
    a = 1 : yes!  set a to 1
  }
}
As always, when adopting anyone else's code, be sure to test it to make sure it works as described.
myreil
Posts: 11
Joined: Mon Mar 14, 2011 2:25 pm

Re: Time dependent conductance in DERIVATIVE block

Post by myreil »

Dear Ted,
Thank you very much for your reply!!

But how should I call it from hoc level? I have something like:
objref nc
soma nc = new NetCon(&v(0.5),target,5,0,0)
where the target is the mod file with the NET_RECEIVE command.
I get the following error:
{run()}
^
finitialize(-65)
init()
stdinit()
run()

Do you know what could possibly cause this error?

Actually what I am trying to do is to have as a source a pulse or a series of pulses. I want to have something like IClamp, but to define it myself. By the time the pulse in on the conductance is changing and according to the kinetic model to get the current.
Would this be possible?

To get this straight , I want to incorporate to an hh model a new ion channel, which is activated by light. So to begin with I thought to have a simple cell and trigger it by a series of pulses, which would correspond to light. I have the kinetic model of the channel and the conductance as said is changing as time passes and is characterized by a different equation when light is on and when light is off. Th source would be the light and the target an hh model with the new ion channel incorporated in it.

Hope that makes sense and that you could give me some help.
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Time dependent conductance in DERIVATIVE block

Post by ted »

You misunderstand the purpose of the NET_RECEIVE block in this mechanism. Let me explain again.

With fixed step integration, it is OK for NMODL code to use conditional statements to change parameter values in the course of a simulation. Comparisons such as
if (t>t0) . . .
always work as expected with fixed dt.

But with adaptive integration, the integrator chooses the time step. The time step can be so big that t0 will be missed--and the resulting error could be quite large. This is where events are useful. You wanted a mechanism in which a variable was defined by one function for t smaller than some t0, and by different function when t was larger than t0 (or whatever you called the time at which the change was to occur). I suggested modifying the INITIAL block so that the mechanism generates a self event at the start of the simulation. The self event returns at time t0 (or whatever you called it). The purpose of the NET_RECEIVE block is to tell the mechanism what to do when the self event arrives. In the example I provided, all that happened was that an ASSIGNED variable called a was given a different value than it started with. (I realize now that the INITIAL block should assign some other value to a, e.g. 0). So based on that toy example, how would you modify FUNCTION Ga2 from the form that works with fixed dt? Hint: in the if (...) statements, instead of comparing t with t0, compare a with 0 or 1.

NetCons have nothing whatsoever to do with such usage of the NET_RECEIVE block.

But now that I think about it, I remember that density mechanisms can't have NET_RECEIVE blocks. Tell me, is your mechanism a density mechanism (has a SUFFIX statement in the NEURON block) or a point process? If it's a density mechanism, I'll have to come up with a different suggestion for you.
myreil
Posts: 11
Joined: Mon Mar 14, 2011 2:25 pm

Re: Time dependent conductance in DERIVATIVE block

Post by myreil »

Dear Ted,

Thank you so much! I figure it out how it works and what my problem is!! My process is a POINT PROCESS so I can put the NET_RECEIVE command.
The problem actually is that my variable Ga1 changes with time and one equation of t characterizes it when ton<t< toff and a different one when t>toff. So in our toy example something like this:
INITIAL {
net_send( 1, toff)
}
NET_RECEIVE (w){
if (flag==1){
Ga2=a*t
} else{
Ga2=b*t+c
}


So in the case of NET_RECEIVE it will detect the event at toff and change the value of Ga2 but how can this continuously change with time??
In a simple case a FUNCTION would work, but as I am using Ga2 in the DERIVATIVE block it will be unstable.

Do you have any idea how could i make something like this work?

I really appreciate your help!!!
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Time dependent conductance in DERIVATIVE block

Post by ted »

myreil wrote:Ga1 changes with time and one equation of t characterizes it when ton<t< toff and a different one when t>toff
What happens when t<ton?
So in our toy example
Your toy example, not mine. Yours won't do what you want because the statements in a NET_RECEIVE block are evaluated only when an event arrives. The point process will receive only one event: the self event that it launched at initialization.

The key to the solution is suggested in my two posts from May 31. Let me be more explicit.

In your conceptual model, Ga2 is described by f1 for t<toff, but by f2 when t>=toff. It is not a STATE, and it doesn't even have to be an ASSIGNED. If you are always going to use fixed dt, you could simply write

Code: Select all

FUNCTION Ga2(t) {
  if (t<toff) {
    Ga2 = f1(t)
  } else {
    Ga2 = f2(t)
  }
}
where f1(t) and f2(t) are algebraic expressions involving t (replace them with whatever you want).

So instead you use a self event that changes the value of a new ASSIGNED parameter (for this example call it a). Initialization sets a to 0 and launches a self event that returns at toff (your mechanism has a DERIVATIVE block with three STATEs O1, O2, and C1 whose values must be initialized by the INITIAL block). Arrival of a self event sets a to 1. The function for Ga2 could then be

Code: Select all

FUNCTION Ga2(t) {
  if (a==0) {
    Ga2 = f1(t)
  } else {
    Ga2 = f2(t)
  }
}
and this should work during a simulation. Your DERIVATIVE block would be

Code: Select all

DERIVATIVE states {
   O1'=Ga1(t)-(Gd1+e_ct+Ga1(t))*O1+(e_ct-Ga1(t))*O2-Ga1(t)*C1
   O2'=e_ct*O1-(e_tc+Gd2)*O2+Ga2*C1
   C1'=Gd2*O2-(Gr+Ga2)*C1
}
But this wouldn't allow you to use a simple for loop like

Code: Select all

for i=0,100 {
  t = i*0.1
  print t, objrefname.Ga2(t)
}
to examine how Ga2 varies with t. Why? Because a is what determines whether the value of Ga2 is generated by f1 or f2, and a will not change unless you run a simulation.

So the function for Ga2 should be

Code: Select all

FUNCTION Ga2(t) {
  if ((a==0) && (t<toff)) {
    Ga2 = f1(t)
  } else {
    Ga2 = f2(t)
  }
}
and this will work during a simulation, or if you call the function directly from hoc with a simple for loop (as long as you call finitialize() first so that a is set to 0), e.g.

Code: Select all

finitialize()
for i=0,100 {
  t = i*0.1
  print t, objrefname.Ga2(t)
}
myreil
Posts: 11
Joined: Mon Mar 14, 2011 2:25 pm

Re: Time dependent conductance in DERIVATIVE block

Post by myreil »

ted wrote:
myreil wrote:Ga1 changes with time and one equation of t characterizes it when ton<t< toff and a different one when t>toff
What happens when t<ton?

When t<ton Ga2=0 and it changes for the first time when t>ton.

Ok, I think I understood everything you have said so far and thank you again very much! However I still get the error:

oc>/Applications/NEURON-7.1/nrn/umac/bin/nrniv.app/Contents/MacOS/nrniv: Bus error See $NEURONHOME/lib/help/oc.help
near line 10
{run()}
^
fadvance()
advance()
step()
continuerun(150)
and others

where 150 is tstop. ( If I set tstop<ton, which doesnt really helps but I just tried it, it works)

if I set ton=0 then i get the following error

oc>/Applications/NEURON-7.1/nrn/umac/bin/nrniv.app/Contents/MacOS/nrniv: Bus error See $NEURONHOME/lib/help/oc.help
near line 10
{run()}
^
finitialize(-65)
init()
stdinit()
run()

What causes this problem and how could I fix it?

Thank you again so much!!
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Time dependent conductance in DERIVATIVE block

Post by ted »

I see that there were several typographical errors in my previous message. Each instance of

Code: Select all

FUNCTION Ga2(t) = {
should have been

Code: Select all

FUNCTION Ga2(t) {
That should have generated an error that would stop compilation of the mod file--so I bet you didn't simply copy my mistake! I have now gone back and fixed the error.
When t<ton Ga2=0
Then, in addition to initializing the STATEs, your INITIAL block also needs to set Ga2 to 0, and the self event that it sends should be used to switch Ga2 from 0 to the form it must have when ton<t<toff.

Code: Select all

INITIAL {
  a = 0
  Ga2 = 0
  net_send(ton, 1)
  . . . initialize STATEs . . .
}
And the NET_RECEIVE block needs to respond to arrival of the { self event with flag 1 } by launching a { self event with flag 2 }.

Code: Select all

NET_RECEIVE () {
  if (flag==1) {
    a = 1
    net_send(toff-ton, 1)
  }
  if (flag==2) {
    a = 2
  }
}
And the FUNCTION has to be changed so that Ga2 is calculated properly.

Code: Select all

FUNCTION Ga2(t) {
  if ((a==0) && (t<ton)) {
    Ga2 = 0
  } else {
    if ((a==1) && (t>=ton) && (t<toff)) {
      Ga2 = expression 1 involving t
    } else {
      Ga2 = expression 2 involving t
    }
  }
}
Regarding the error message you encountered, I suggest you update to a more recent version of NEURON. For OS X, that means either getting the most recent development source code from the mercurial repository and compiling it, or trying the most "recent" version 7.2 dmg file, which I see dates all the way back to June 10 last year--
http://www.neuron.yale.edu/ftp/neuron/v ... 0.4.11.dmg
Post Reply