Myelinated Nerve with Extracellular Potential

[Matte]

Hey all,
first of all, thanks a lot for the great forum, I have been reading for quite some time now, and got pretty far with my model, and just have a few question left. My model is basically the model of the myelinated nerve by McIntyre 2002 () for which I added an extracellular potential to excite the nerve. I am using the python interface to set up all the parameters of the fiber, node number, electrical properties and so on, and implemented a matlab-python interface to be able to call the Neuron solver from Matlab, nothing really fancy...
The model works right now, but I just wanted to make sure that I have implemented everything correctly, because I am relatively new to the field of modelling nerve fibers.

1. Setting up the extracelluar potential:
After having declared a couple of parameters, I am basically performing the following lines of codes to create and connect the different segments of the nerve fiber. This is done in python using the h("statement") approach.

Code: Select all

for i=0,axonnodes-1 {
		node[i]{					
			nseg=1
			diam=nodeD
			L=nodelength
			Ra=rhoa/10000
			cm=2
			insert axnode			
			insert extracellular xraxial=Rpn0 xg=1e10 xc=0
			}
		}
	for i=0, paranodes1-1 {
		MYSA[i]{
			nseg=1
			diam=fiberD
			L=paralength1
			Ra=rhoa*(1/(paraD1/fiberD)^2)/10000
			cm=2*paraD1/fiberD
			insert pas
			g_pas=0.001*paraD1/fiberD		
			e_pas=-80
			insert extracellular xraxial=Rpn1 xg=mygm/(nl*2) xc=mycm/(nl*2)
			}
	}
	for i=0, paranodes2-1 {
		FLUT[i]{
			nseg=1
			diam=fiberD
			L=paralength2
			Ra=rhoa*(1/(paraD2/fiberD)^2)/10000
			cm=2*paraD2/fiberD
			insert pas
			g_pas=0.0001*paraD2/fiberD		
			e_pas=-80
			insert extracellular xraxial=Rpn2 xg=mygm/(nl*2) xc=mycm/(nl*2)
			}
	}
	for i=0, axoninter-1 {
		STIN[i]{
			nseg=1
			diam=fiberD
			L=interlength
			Ra=rhoa*(1/(axonD/fiberD)^2)/10000
			cm=2*axonD/fiberD
			insert pas
			g_pas=0.0001*axonD/fiberD
			e_pas=-80
			insert extracellular xraxial=Rpx xg=mygm/(nl*2) xc=mycm/(nl*2)
			}
	}
	for i=0, axonnodes-2 {
		connect MYSA[2*i](0), node[i](1)
		connect FLUT[2*i](0), MYSA[2*i](1)
		connect STIN[6*i](0), FLUT[2*i](1)
		connect STIN[6*i+1](0), STIN[6*i](1)
		connect STIN[6*i+2](0), STIN[6*i+1](1)
		connect STIN[6*i+3](0), STIN[6*i+2](1)
		connect STIN[6*i+4](0), STIN[6*i+3](1)	
		connect STIN[6*i+5](0), STIN[6*i+4](1)	
		connect FLUT[2*i+1](0), STIN[6*i+5](1)
		connect MYSA[2*i+1](0), FLUT[2*i+1](1)
		connect node[i+1](0), MYSA[2*i+1](1)	
		}

In the first three loops, this codes also adds extracelluar mechanisms. This was also done in the code provided by the author, but it is especially required to be able to apply a stimulus via the extracellular potential (please correct me, if I am wrong). I am then defining the extracelluar potential at all nodes (i.e., nodes of Ranvier and myelin segments MYSA, STIN, FLUT) and modulate this by a stimulus waveform (such as a sin wave for the simplest case).

Code: Select all

# create vector of time-varying extracellular potential
v1 = h.Vector( data_A * data_V[1] )
# play into first node of Ranvier
v1.play(h.node[1](0.5)._ref_e_extracellular, h.dt)

data_A refers to the modulation waveform, data_V[1] to the extracellular potential. This is done for all nodes and myelin segments.
My first question:
Is this correctly implemented? Do I have to change or add the default values following insert extracellular (i.e., xraxial, xg, xc) in order for the extracellular potential to function correctly?

2. Boundary Handling:
As far as I know, the default handling of the ends of the fiber is to use sealed ends (no axial current flow). I have experienced that this results in the fact that action potentials are mostly initiated at the boundary nodes. What I would actually like to implement is something like "solely axial current flow", so basically I'd like to model a finite part of an infinitely long fiber.

My second question:
Is there a way to implement such a behaviour? I checked the documentation/forum, but I didn't really find an answer to this problem.

Thanks a lot in advance for any helpful comments! :-)
Warm regards,
Matte

[ted]

There's an old saying that goes something like "the only thing than debugging your own code is debugging somebody else's code." There's another old saying to the effect that "the only thing more treacherous than a program that crashes is a program that works" (that is, a program with a bug that does not force a failure). And finally there's the general advice against reinventing the wheel. Which adds up to this: even though I may not see an "error" as such in the code you present, and your program may not generate a syntax or a runtime error when fed to NEURON, that doesn't guarantee that your program is "correct".

There are two primary concerns. First, is the model implemented by your Python code identical to the hoc original. The first concern could have been largely avoided by simply reusing hoc code lifted straight from the implementation by McIntyre et al.. After all, anything in a hoc-specified model can be accessed from Python. But either reusing McIntyre's code, or re-implementing the model with Python, still leaves you with the necessity of verifying that what you have actually replicates the original model. That can be done by using NEURON's Model View tool to analyze both the original and your new implementation, to verify that the anatomical and biophysical parameters are identical, and also checking that the topology is identical (the latter can be done by comparing the output generated by executing topology() or h.topology(), depending on whether the interpreter is hoc or Python). Another important test is to verify that the new implementation generates numerical results that are identical to those produced by the original.

The second concern is whether your program implements extracellular stimulation correctly, that is, in a way that produces an extracellular potential that varies along the axon in the way you expect it to vary. Only you can discover this. Plot it as a function of distance along the axon and see if it matches your expectation.

Do I have to change or add the default values following insert extracellular (i.e., xraxial, xg, xc) in order for the extracellular potential to function correctly?

No, the defaults effectively isolate each node just outside the membrane of each segment from the nodes of the adjacent segments, while providing a direct connection between that node and its corresponsing e_extracellular voltage source.

What I would actually like to implement is something like "solely axial current flow", so basically I'd like to model a finite part of an infinitely long fiber.

Just make the model axon long enough that the area of interest is far from the region(s) in which sealed end effects are significant. How long? Plot peak spike amplitude vs. distance along the model. How big is a "significant" effect? That's your decision--do something that you can explain to those who may review whatever publications will result from your work.