Applying Sinusoidal Voltage to Membrane

[dchu713]

Hello!

I am trying to apply a sinusoidal voltage to vary the membrane potential over time of a single axon model and I am would like some help figuring out the best way to do it. I am considering a VClamp in hoc, but I also saw the sinusoidal IClamp (viewtopic.php?p=19872&hilit=SEClamp#p19872) and have been looking through how to apply that. I have also been exploring the extracellular time varying input thread (viewtopic.php?p=8318&hilit=pulse+amplitude#p8318). I am very much a novice in this field, so I may have missed the easiest application!

I am using the Gaines et al (2018) paper (https://modeldb.science/243841?tab=2) as inspiration and for the .hoc and mod files. Based on my understanding of the codes, the extracellular values of every nodal, paranodal, and internodal segment is initialized as zero until the time equals the delay. During the length of time equal to the pulse width, the extracellular values are set to the values that were determined by ANSYS MAXWELL. I used a COMSOL multiphysics model to find the extracellular values and changed the .dat files to what I calculated. The overall syntax of the .hoc file is largely the same.

I tried to loop the extracellular values for the duration of the pulse width, but with my new understanding of the code, I do not think that is the correct way to do it. I also tried using the vector.sin function to generate a sinusoidal function within the extracellular voltage vectors, but I believe this varied the values over space and not time. Below is the part of the code that I am working with the most.

Code: Select all

//the following procedure only gets called at t=0, t=delay, t=delay+pulsewidth//
// changes extracellular stimulation values at the beginning of of the test and at the beginning and end of the stimulation pulse
proc set_e_ext() {

    // beginning of the test
	// extracellular voltage is zero, wait for t = delay to start pulse

	/* 
		This first conditional sets all the extracellular values to 0 and waits until the simulation time reaches the defined
		delay time. I have the "no pulse" print to see what stage its at when running
	*/  
	if ($1 == 0) {
		forall e_extracellular = 0//set extracellular voltage source to ground
		cvode.event(delay, "set_e_ext(1)")  // when to turn on
		print "no pulse"

    // during the stimulation pulse
	// set extracellular voltage to values determined by COMSOL
	// wait for the end of the pulse (t = delay + pw)

	}else if ($1 == 1) {
		print "pulse now!"
		for t=delay,delay+pw { // I added this loop condition to see what would happen ~ dchu713
		for i=0,axonnodes-1 {
			node[i]{
			e_extracellular = node_volt.x[i]			
				}
		}
		for i=0, paranodes1-1 {
			MYSA[i]{	
			e_extracellular = mysa_volt.x[i]
				}
		}
		for i=0, paranodes2-1 {
			FLUT[i]{
			e_extracellular = flut_volt.x[i]
				}
		}
		for i=0, axoninter-1 {
			STIN[i]{
			e_extracellular = stin_volt.x[i]	
				}
		}	
		}	
//This cvode.even resets the extracellular voltages to zero after the pulse width
		 cvode.event(delay + pw, "set_e_ext(2)")
    // end of the stimulation pulse
	// reset extracellular stimulation to zero
	// wait for the beginning of a new test

	}else if ($1 == 2) {
		print "Other case?"
		forall e_extracellular = 0
		cvode.event(0, "set_e_ext(0)")
	}
}

Thanks for any support!

[ted]

First you need to decide exactly what you are trying to do.

This
"apply a sinusoidal voltage to vary the membrane potential over time"
means you want to force membrane potential to follow a specific time course. That is best done with a voltage clamp--specifically, an instance of NEURON's SEClamp class driven by a sinusoidal forcing function. VClamp would not be the tool to use; VClamp was implemented specifically to study the effects of imperfect instrumentation on (squid axon) voltage clamp experiments. And the "sinusoidal IClamp" discussed in the post that you cite would only inject a sinusoidal current, and that would not make membrane potential follow a sinusoidal time course (because of channel gating).

But an SEClamp injects current into a single segment (compartment), so this
"of a single axon model"
means that each segment of your model axon would need its own instance of the SEClamp class, and all of these SEClamps must be driven by some forcing function.

Then in the next paragraph it looks like you might actually want to force extracellular potential to follow a particular time course. And that would not be done with voltage clamps or current clamps.

So what is it that you want to do? Force membrane potential to follow a particular time course, or force extracellular potential to follow a particular time course?

[dchu713]

Hello Ted,

Thank you for your reply! After some review, I am trying to force the extracellular potential to follow a particular time course. After all, I can always try both methods and see the results! Since I am using the Gaines 2018 paper as a reference and they used the extracellular voltage values, I would like to try and force my extracellular values to follow the time course. That way I can best compare the differences between the values.

Thank you so much for your questions, they really helped!

Best,
dchu713

[ted]

OK, three more questions, and then we can move on to specific suggestions for you.

Are you assuming that the extracellular medium is purely resistive?
Do you already have the time course of extracellular potential that you want to apply to the axon?
Which programming language do you prefer to use in your project: hoc or Python?

[dchu713]

Hello Ted,

Thank you for your continued patience and support! I have the answers to your questions

Are you assuming that the extracellular medium is purely resistive?

Yes, we can assume that the extracellular medium is purely resistive.

Do you already have the time course of extracellular potential that you want to apply to the axon?

I do have a time course in mind. I want to apply a sine pulse with a amplitude of 1 V. I am trying to apply the voltage from an anode-cathode pair in the center of the model and have that pulse propagate down the length of the model. I have extracellular voltages calculated at the center point of each segment using a COMSOL model in a similar manner to Gaines et al (2018). Essentially what I am trying to do is have each value vary with a pulse over time.

Which programming language do you prefer to use in your project: hoc or Python?

I prefer hoc for this project, though I am not opposed to learning Python if it is more applicable or if it is a better way to do this.

I hope that clarifies the problem! If not, please ask any necessary questions.

Thank you again and best regards,
dchu713

[ted]

I am trying to force the extracellular potential to follow a particular time course. . . .
I want to apply a sine pulse with a amplitude of 1 V.

So the electrode-tissue interface is not a part of your conceptual model. Instead, it is assumed that it is possible to force extracellular potential to a particular value.

I have extracellular voltages calculated at the center point of each segment using a COMSOL model

And the question is how to use these values to drive your computational model.

. . . we can assume that the extracellular medium is purely resistive. . . .

and is also presumably linear. This has two important consequences.

(1) Linearity implies that the extracellular potential at any point is directly proportional to the current supplied by the stimulator. Furthermore, in a "purely resistive" conductive medium, the relationship between stimulus current and extracellular potential is instantaneous. That is, given a stimulus current with time course i(t), and some location p, then vext(p, t) (i.e. the extracellular potential at p as a function of time) is given by
vext(p, t) = r(p) * i(t)
where r(p) has units of resistance and has sign and magnitude that vary with the location of p.

(2) From (1) it follows that the extracellular potential waveform at any point p2 is the same as that at any other point p1 multiplied by a constant whose magnitude and sign depend on the locations p2 and p1. After all, the time courses of vext at p2 and p1 are
vext(p2, t) = r(p2) * i(t)
and
vext(p1, t) = r(p1) * i(t)
so
vext(p2, t)/vext(p1, t) = r(p2)/r(p1)
vext(p2, t) = k12 * vext(p1, t)
where
k12 = r(p2)/r(p1)

Note that k12 depends only on the locations of p1 and p2 in space.

With this background, you are ready to see how to use a density mechanism called xstim which I developed in order to implement extracellular stimulation in a resistive medium. https://modeldb.science/239006 contains source code that shows how I used it to model extracellular stimulation of an axon. You should also read the paper that is associated with this model (let me know if you have difficulty finding that paper).

xstim has a POINTER variable called ex that is used to drive e_extracellular. The hoc code that couples each segment's ex_xstim to its e_extracellular is

Code: Select all

forall {
  insert xstim
  for (x,0) setpointer ex_xstim(x), e_extracellular(x)
}

(the iterator for (x,0) omits the zero-area nodes at 0 and 1).

This block of code inside xstim.mod

Code: Select all

BEFORE BREAKPOINT { : before each cy' = f(y,t) setup
  ex = is*rx*(1e6)
}

ensures that the values of ex throughout the model are calculated by NEURON's computational engine before the BREAKPOINT blocks are executed.

is is a GLOBAL ASSIGNED variable in xstim, but xstim has no statement that assigns a value to is. Instead, is must be driven by something outside of xstim. That "something" could be a Vector (i.e. by Vector.play), but I used a couple of NMODL-specified function generators: Fsquare (a rectangular pulse generator) or Fzap (a swept sine wave generator). These function generators have a POINTER variable x that is coupled to is_xstim by a statement like

Code: Select all

setpointer fzap.x, is_xstim

or

Code: Select all

setpointer fsq.x, is_xstim

Note that is_xstim does not reference a particular location in a particular section, and the setpointer statement is not contained in an iterator that loops over all segments in all sections. Instead, since is_xstim is a GLOBAL variable, only a single setpointer statement is needed, and that statement will affect all segments in all sections that have xstim.

"This is all great, but I'm not dealing with stimulus current. I'm dealing directly with extracellular potential."

Of course. Note that ex_xtra must be in mV. In xstim.mod, ex is calculated by
ex = is*rx*(1e6)
where the units of ex, is, and rx are mV, mA, and megohms. For your particular application I recommend that you set the amplitude of is to 1 (do that by setting the amp parameter of Fzap or Fsquare to 1). Then you can use the value of rx, which is a RANGE variable, to ensure that each segment's ex has the correct numerical value.

How?

If the numerical value of rs is set equal to the maximum local extracellular potential in V, then the product
rx*(1e6)
will be the numerical value of the maximum local extracellular potential in uV--a number that is 1000 times too large. In order to get the numerical value of the local peak voltage in mV, you need to set
rx = 1e-3 * maximum local extracellular potential in V

So, in pseudocode, here's what you need to do:

Code: Select all

for each section in the model
  for each segment in this section
    set rx = 1e-3 * the COMSOL-calculated maximum local extracellular potential
    (in V) at this particular location

And since you're using a bipolar stimulus, the "maximum local extracellular potential" should have the correct sign. For those locations that go positive when the anode is delivering positive current, the sign should be positive; for locations that go negative when the anode is positive, the sign should be negative. The easy way to discover sign is to use COMSOL to apply a single "positive stimulus pulse" and see which locations go positive and which go negative.

[dchu713]

Hi ted,

Thank you so much for all the info! I am able to access the paper that you linked so I will read through that and I'll implement xstim into my model and that should drive the extracellular voltage to follow my time course.

Thank you again!

Best,
dchu713