Extracellular Stimulation (2d)

[AndiL]

Hey Guys!

Im writing on my thesis and need a chapter about extracellular stimulation in neuron

i read some papers about that topic, all the forum-threads and looked into the examples provided (extracellular_stim_and_rec and ecstim).
the problem is that i think i just dont get it right - still...

lets say i need (as a typical scenario)
21 nodes combined with 20 internodes (or just 2 nodes with 1 internode for a start, or even just ONE node only)
an electrode near the center node (or whereever in the x/y space)
all that in a homogenous liquid (simulated via extracellular)

attempt one:
it should be enough to re-design the "cell.ses" which i can do with the topology manager...
this would be a solution - but i still think that a 2d solution might be the better attempt...


im working and trying now since a while and end up a little frustated.
can anyone in here help me, please?

thankful for every suggestion.
cheers,
andi

[ted]

By "2d" do you mean "second" (as in "alternative") or do you mean "two dimensional"?

[ChrisR]

Hello everybody!

Actually this is the simple case I've started with before messing around with the suction electrode configuration. Maybe this is the configuration you are looking for.

So this should work for the myelinated frog axon configuration using the files interpxyz.hoc, setpointers.hoc, calcrxc.hoc and stim.hoc povided by http://www.neuron.yale.edu/ftp/ted/neur ... nd_rec.zip. Just be sure to set your extracellular rho (resistivity of tissue) in calcrxc.hoc (e.g. 300)

But I am not totally sure if the extracellular configuration is correct:

for nodes: insert extracellular xraxial=1e9 xg=1e9 xc=0 e_extracellular=0
for internodal sections: insert extracellular xraxial=1e9 xg=0 xc=0 e_extracellular=0

Ted, what do you think?

Code: Select all

load_file("nrngui.hoc")

tstop=5
v_init = -70
celsius = 14
d = 5 // axon diameter

rho_i = 110 // Ohm cm (Stampfli et al. 1952)
rho_e = 300 // Ohm cm (Abzug et al. 1974)

cm_n = 2 // uF/cm2 (Frankenhaeuser and Huxley 1964)
gm_n = 0.0303 // S/cm2 (Frankenhaeuser and Huxley 1964)

cm_m = 0.005 // uF/cm2 (Moore JW, Joyner RW, Brill MH, Waxman SD, Najar-Joa M (1978) Biophys J 21:147-60)
gm_m = 1.5e-6 // 1.5e-6 S/cm2 (Moore JW, Joyner RW, Brill MH, Waxman SD, Najar-Joa M (1978) Biophys J 21:147-60)

l = 2.5 // um (Dodge and Frankenhaeuser 1959)
rifd = 100 // ratio of internodal space to fiber diameter (Dodge and Frankenhaeuser 1959)
rafd = 1/1.5 // ratio of axon to fiber diameters (Goldman and Albus 1968) (Dunn 1977, cited by Honrubia 1989)

Df = d/rafd 
Li = Df*rifd // length of internodal section 

nn=21 // number of nodes
Er= -70 // mV

create myelin[nn-1], node[nn]

for i=0, nn-2 {
    connect myelin[i](0), node[i](1)
    connect node[i+1](0), myelin[i](1)
}

forsec "node" { 
    L = l
    diam = d
    Ra = rho_i    
    cm = cm_n
    
    insert fh
    nai=13.74  nao=114.5  ki=120 ko=2.5
    
    insert extracellular xraxial=1e9 xg=1e9 xc=0 e_extracellular=0
    insert xtra 
}

forsec "myelin" {
    L = Li
    diam = d
    Ra = rho_i 
    cm = cm_m
   
    insert pas
    e_pas = Er
    g_pas = gm_m 

    insert extracellular xraxial=1e9 xg=0 xc=0 e_extracellular=0
}

forsec "myelin" { nseg = 1 }
forsec "node" { nseg = 10 ek = -90}

forall {
    define_shape()
}

load_file("interpxyz.hoc")      // only interpolates sections that have extracellular
load_file("setpointers.hoc")    // automatically calls grindaway() in interpxyz.hoc
load_file("calcrxc.hoc")        // computes transfer r between segments and recording electrodes
load_file("stim.hoc")          // extracellular stimulus

access node[10]
xpos = node[10].x_xtra(0.5) // position of center node
setelec(xpos, 1000, 0)

[ted]

I should have been more explicit about usage in writing the readme.txt that accompanies the files in extracellular_stim_and_rec.zip. Of course, at the time it was all very clear to me, which is why some helpful tips were left unstated. Here I try to atone for these omissions.

First, the code is written in a very modular fashion, which separates specification of biological properties, instrumentation, control code, and user interface from each other. This was done for several reasons, not least of which is the fact that it allows "plug and play" substitution of one model cell for another.

To see what this means, consider the file initxstim.hoc:

Code: Select all

load_file("nrngui.hoc")
load_file("cell.ses")
load_file("anatscale.hoc")      // show xyz scale bars
load_file("interpxyz.hoc")      // only interpolates sections that have extracellular
load_file("setpointers.hoc")    // automatically calls grindaway() in interpxyz.hoc
load_file("rigc.ses")           // RunControl, graphs of v, vext, e_extracellular
load_file("rigxc.ses")          // graph of er_xtra, just for diagnostic|development purposes
load_file("field.hoc")          // computes extracellularly recorded potential vrec
load_file("vrecc.ses")          // graph of vrec(t)
load_file("calcrxc.hoc")        // computes transfer r between segments and recording electrodes
// load_file("iclampc.ses")     // intracellular stimulus
load_file("stim.hoc")           // extracellular stimulus
load_file("msgxc.hoc")          // minimal hints

The second load_file() statement tells NEURON to read the file that specifies the anatomical and biophysical properties of the model cell. In this particular case, the file is cell.ses, which was created by saving a CellBuilder to a ses file, but it could just as easily be a hoc file. The file must do the following--
1. Specify topology: create and connect the sections, and use a single access statement to declare one of them as the default section.
2. Assign biophysical properties and specify discretization: assign values to Ra, cm and nseg; insert whatever density ("distributed") mechanisms are necessary.
3. Insert extracellular and xtra into all sections that are to "feel" the effect of an extracellular field. The parameters of extracellular should be assigned the following values:
xraxial=1e+09 xg=1e+09 xc=0 e_extracellular=0
4. If geometry was specified with the stylized method (i.e. L and diam) rather than the pt3d method, call define_shape(). Once is enough.

The two other files that may require changes by the user are calcrxc.hoc and stim.hoc.

The introductory comments to calcrxc.hoc are somewhat misleading. The code does not set up bipolar stimulation--it sets up monopolar stimulation. I offer my deepest apologies for this remnant of a previous implementation and promise never to do such a thing again. Clever users will quickly discover how to make it generate a bipolar stimulus (hint: the solution is in the file itself).

stim.hoc contains the code that specifies the waveform of the current that is applied to the stimulating electrodes.

<boilerplate>
Conditions of reuse

Anyone is welcome to reuse any or all of it and modify it as necessary, but if such reuse produces fame or fortune I would appreciate a piece of the action. It would be nice if results published in a reputable scientific journal were accompanied by acknowledgment of the origins of the code. And of course I wrap myself in the same litany of disclaimers that MicroSoft repeats in its EULAs.

I'd be glad to hear of suggestions for improving this code. Support is conditioned on the requirement that the user's code adhere to good programming practice. ("well, then why are we using hoc?")
</boilerplate>


ChrisR, you did a good job of deciphering the basic structure of initxstim.hoc and identifying the files that are essential for extracellular stimulation.
You will want to change the parameters of extracellular as noted above.
xtra should be inserted into all sections of your model, unless the myelinated internodes are encased in little Faraday cages.

Here are some other questions and comments that may help you avoid some nasty little gotchas:
Is there any particular reason for the nseg values you're using, i.e. 1 for the (long) internodes, and 10 for the (very short) nodes? Also, when nseg is even, what does it mean to refer to somerangevariable(0.5) (or even somerangevariable without the (0.5) bit)? It would be a good idea to read the FAQ list's entries about nseg -- see http://www.neuron.yale.edu/neuron/faq/general-questions.
Finally, please see 'Use only one "access" statement' in the Forum's "Hot tips" area viewforum.php?f=28

[AndiL]

Im using the Code as ChrisR posted it... i think it works out for me (changed the models from McNeil in ModelA, B and C described by Richardson as predecessor of the MRG Modell)

However i have a couple of questions:
1) if i want to change the eletrical properties of the substance in which is surrounding the fiber it should be enough to change the values of extracellular, right?

2) why do i have to put a negative amplitude of the Stimulus?

3) if i try the code with the MRG Modell i get a AP but a kind of outburst of endless APs following. Any idea why this is happening?

well id appreciate any comments and answeres...
thanks in advance,
andi

[ted]

Sorry, but I'll need some more information to answer your questions.

1) if i want to change the eletrical properties of the substance in which is surrounding the fiber it should be enough to change the values of extracellular, right?

What electrical properties of the conductive medium do you want to change, and which "values of extracellular" do you propose to change?

2) why do i have to put a negative amplitude of the Stimulus?

What are you referring to? Can you point to a specific line or block of source code?

3) if i try the code with the MRG Modell i get a AP but a kind of outburst of endless APs following. Any idea why this is happening?

What's the MRG model?

[AndiL]

1) if i want to change the eletrical properties of the substance in which is surrounding the fiber it should be enough to change the values of extracellular, right?

What electrical properties of the conductive medium do you want to change, and which "values of extracellular" do you propose to change?

well i want to adopt the parameters of the medium that they fit the one of lets say: fat (e.g.) cant i do that by changing the xg in the insert extracellular part?

2) why do i have to put a negative amplitude of the Stimulus?

What are you referring to? Can you point to a specific line or block of source code?

if i run the program, i have box called "extracellular stimulus current" where i can define the parameters for the electrode (such as delay, duration, amplitude)
standard here is del=1, dur=2, amp=-0.3
if i change the amp to smaller values (-0.5, -1.0, ...) it works out as it should.
if i change it to positive values it kind of goes... crazy. meaning: showing values extremly negative. (like -40.000mV)

3) if i try the code with the MRG Modell i get a AP but a kind of outburst of endless APs following. Any idea why this is happening?

What's the MRG model?

MRG = McIntyre, Richardson, Grill Modell (or just McIntyre) which is already a double cable model (using extracellular for modelling the myelin sheath)
first i thought that maybe the double use of extracellular (as myelin AND medium) might be the problem.
HOWEVER: ModelC does the same. and the only difference here is the caluclation of the ion-flows (dont think that this is the problem) AND the seperation of the internodal compartments in 3 different parts (MYSA, FLUT and STIN).
so could it be that because of the seperation of this internodal compartments if creates the outburst?
so what i thought is: is there a possability to but the extracellular (medium) over ALL 3 of the internodal comparments at once? (and not over every single one?)

thanks for every comment or thought spared

[ted]

i want to adopt the parameters of the medium that they fit the one of lets say: fat (e.g.) cant i do that by changing the xg in the insert extracellular part?

Only if the geometry of the extracellular medium is more or less cylindrical and can be represented by a few concentric layers.

2) why do i have to put a negative amplitude of the Stimulus?
. . .
if i run the program, i have box called "extracellular stimulus current" where i can define the parameters for the electrode (such as delay, duration, amplitude)
standard here is del=1, dur=2, amp=-0.3
if i change the amp to smaller values (-0.5, -1.0, ...) it works out as it should.
if i change it to positive values it kind of goes... crazy. meaning: showing values extremly negative. (like -40.000mV)

First, a comment: I have never seen or heard of anything like this happening.

Now, some questions to try to figure out what's going on.

-40 mV is extremely negative? Do you mean -4*10^4 mV?

Did you place the stimulating electrode right on top of, or very close to, a section? The potential generated by a point source of current in an infinitely large conductive medium is proportional to 1/r^2, where r is the distance from the point. If r is very small, you should see a large potential.

3) if i try the code with the MRG Modell i get a AP but a kind of outburst of endless APs following. Any idea why this is happening?

Nope. I don't see anything like that myself when I take their biological specification and instrument it with my extracellular stimulation code.

Which is not entirely a trivial undertaking.
1. MRGaxon.hoc contains the statement
load_proc("nrnmainmenu")
which is best to omit for reasons noted in the Programmer's Reference's comments about load_proc.
2. The intracellular stimulus specified at the end of MRGaxon.hoc should be commented out.
3. load_file("cell.ses") must be replaced by
load_file("MRGaxon.hoc")
forall insert xtra
define_shape()
(the latter is essential because MRGaxon.hoc uses the stylized (L, diam) specification of geometry, but the strategy employed by initxstim.hoc requires the pt3d data to exist).
4. rigc.ses and rigxc.ses need to be rebuilt because their graphs involve paths with sections that no longer exist.
5. The initial location of the stimulating electrode is (50, 0, 0), which lies right on top of one of the model's sections. As noted above, this will produce unrealistically large local extracellular potentials.
6. The model lies on top of the xyz "geometry cal bars". This is a merely cosmetic problem that has no effect on model properties or simulation results.

I'll zip up my implementation and send it to you.

ModelC does the same.

What's ModelC?

[AndiL]

hey ted!

ModelC does the same.

What's ModelC?

ModelC is a kind of predecessor of MRG (McIntyre). it already worked with the double cable as well but didnt divide the internodal compartments as McIntyre did. so it didnt really match the afterpotentials - but had a way better computation-time.


thanks a lot - the zip you sent me worked out:
i tried a kind of strength-duration threshold testing (on different distances) and it fits with the data published in other papers... =)

now i saw the "bipolar" electrode already mentioned in the calcrxc.hoc file. that could be interesting for me as well and i wanted to ask whether you tried and/or tested that already??
actually what i want to do is: to set an AC on the bipolar electrode(s). did you ever try a set like that?
i mean if the bipolar set in the calcrxc.hoc works - the only thing i need to do is to set the STIM in an AC one, right?

thanks in advance.
greets,
andi

[ted]

thanks a lot - the zip you sent me worked out:

Excellent.

Your next task, if you choose to accept it, is to compare what I sent with the code that didn't work, and report what caused the problems. I would be very interested in knowing whether or not there is an error (introduced by the publisher) in the digital version of The NEURON Book.

i tried a kind of strength-duration threshold testing (on different distances) and it fits with the data published in other papers

That's reassuring.

now i saw the "bipolar" electrode already mentioned in the calcrxc.hoc file. that could be interesting for me as well and i wanted to ask whether you tried and/or tested that already??

It works.

actually what i want to do is: to set an AC on the bipolar electrode(s). did you ever try a set like that?
i mean if the bipolar set in the calcrxc.hoc works - the only thing i need to do is to set the STIM in an AC one, right?

Yep, that's all you have to do.

[AndiL]

actually what i want to do is: to set an AC on the bipolar electrode(s). did you ever try a set like that?i mean if the bipolar set in the calcrxc.hoc works - the only thing i need to do is to set the STIM in an AC one, right?

Yep, that's all you have to do.

is there a elegant why to generate this sinus-vector?
i mean to generate e.g. a "amp" / "time" pair to use in combination with the vectors play method?

stim_amp.sinus(freq, base, dt)
sets the vector to my sinus function (which is already half of all i need)
by multiplying those values i can vary the amplitudes.
i dont fully understand the "dt"-parameter though...

can i get the accurate time-vector fitting to the sinus-vector?
(because for the play mathode i need this one as well...)

so that i really can use the line as implemented in stim.hoc:
stim_amp.play(&is_xtra, stim_time, 1)

[ted]

For this, better than Vector play is to use a mod file that computes the waveform analytically. This works for any frequency and, since the function is computed for the actual time, also works for any dt--even with adaptive integration, in which dt may change in the course of simulation. Later today I'll send you a mod file that produces a sine wave with specified onset time, number of cycles, and amplitude, plus instructions for how to link its output to xtra's "stimulus current" variable.

[AndiL]

i already tried to but the sinus-curve in the xtra.mod File.
I changed the line:
ex = is*rx*(1e6)
to:
ex = is*rx*(1e6)*sin(2*3.141*t*20000

so just adding the "*sin()" part.
20.000 is the frequency (and 2*PI*t for really creating a sinus)
the output shows a sinus-function. (so far so good)

But i wanted to compare it with some with the outcome of a paper (bhadra 2007) i am studying and its not matching...
in this they have 51nodes and put the eletrode (producing the sinus with 10kHz) right over the center node at the distance of 1mm. rho of the squid material is changed to 500 (as specified)
then the put an ICLAMP at the very first node to generate a AP there after 40ms
Then they look at the last node if or not a AP accures (depending on the amplitude of the extracellular-sinus-eletrode)

i tried to reproduce that set, but its just not delivering the same results. i thought that maybe multiplying with the sinus in the xtra.mod doesnt really do the same as working with a "real" sinus-stimulus? (or maybe because its not interpolated?)

hmm... maybe you give it a thought (i can also send you the set i created and just discribed if you would have time to have a look on it it would be great... just tell me).

thank you very much anyway
greets,
andi

[ted]

I changed the line:
ex = is*rx*(1e6)
to:
ex = is*rx*(1e6)*sin(2*3.141*t*20000

Don't forget the closing parenthesis.
Use PI instead of 3.141. NMODL knows that PI is 3.1........
ex = is*rx*(1e6)*sin(2*PI*t*20000)

This starts the oscillating field at t = 0. That may not be appropriate if you are dealing with a model that requires time to settle to a steady state before the field begins to fluctuate.

But i wanted to compare it with some with the outcome of a paper (bhadra 2007) i am studying and its not matching...
in this they have 51nodes and put the eletrode (producing the sinus with 10kHz) right over the center node at the distance of 1mm. rho of the squid material is changed to 500 (as specified)
then the put an ICLAMP at the very first node to generate a AP there after 40ms
Then they look at the last node if or not a AP accures (depending on the amplitude of the extracellular-sinus-eletrode).

Are you quite certain that you have replicated their model, stimulus, and simulation parameters?

You are likely to get quantitatively different results if the specification of your model, or stimulus, or simulation parameters (dt or nseg for example) differ in any way from what they used. Papers generally do not provide sufficient detail to allow exact replication of results. Usually the only way to know all details is to get the original source code from the authors, and even that may not be sufficient--it is not uncommon for authors to temporarily change one or more parameters for a particular figure, and then forget to include that detail in either their paper or their code. Often the best one can achieve is qualitative similarity.

Are Bhadra et al. 2007 using the MRG model?

i tried to reproduce that set, but its just not delivering the same results. i thought that maybe multiplying with the sinus in the xtra.mod doesnt really do the same as working with a "real" sinus-stimulus?

This speculation is false.

(or maybe because its not interpolated?)

This speculation is also false. Interpolation is irrelevant since interpolation is not involved--the mechanism is calculating the local potentials exactly at each time step.

[ted]

Table 2 of Bhadra et al. 2007 lists parameters. I find the following differences between those values and the values used by the MRG axon model that you and I have been working with:
Internodal length 175.2*6 vs. 175.667*6
specific membrane capacitance 2 vs. multiple values (0.66, 1.38, or 2 uf/cm2 in various kinds of sections)
specific myelin capacitance 0.1 vs. 0.00041666... uf/cm2
specific myelin conductance 0.001 vs. 4.166...e-6 S/cm2
axoplasmic resistivity 70 vs. multiple values (147.02794, 642.79155 ohm cm in various kinds of sections)

The internodal length differences are too small to matter. The other differences are quite large. You should expect to see only qualitative similarities to the results published in Bhadra et al. 2007.