setpointer in Python

When Python is the interpreter, what is a good
design for the interface to the basic NEURON
concepts.

Moderator: hines

tfoutz99
Posts: 12
Joined: Wed Dec 12, 2007 3:15 pm

setpointer in Python

Post by tfoutz99 »

Hi everybody,
I am having a hard time implementing setpointers in python.
In hoc, I'm using:

Code: Select all

forall {
	for (x, 0) {
		setpointer im_xtra(x), i_membrane(x)
	}
}
but I can't get h.setpointer to work. Here's what I tried:

Code: Select all

for sec_id in h.allsec():
	sec_id.push()
	for seg in sec_id:
		h.setpointer(seg.im_xtra, seg.i_membrane)
	h.pop_section()
I just get this reply:
TypeError: no HocObject interface for setpointer (hoc type 315)

What am I missing?
-Tom
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: setpointer in Python

Post by hines »

Andrew Davison pointed out to me some months ago that there was no python interface to
setpointer and I sadly forgot about it after considering that it should probably look like
h.setpointer(_ref_hocvar, 'POINTER_name', point_process_object) and
h.setpointer(_ref_hocvar, 'POINTER_name', x, sec=section)
for point and density mechanisms respectively (the latter POINTER_name must have the mechanism SUFFIX appended).

I'll fix this before the end of next week and commit the change to the mercurial repository.
tfoutz99
Posts: 12
Joined: Wed Dec 12, 2007 3:15 pm

Re: setpointer in Python

Post by tfoutz99 »

Thanks!
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: setpointer in Python

Post by hines »

setpointer now works from Python as a function. See:
http://www.neuron.yale.edu/hg/neuron/nr ... fedb04513a

The syntax for setting a POINTER in an NMODL density mechanism is different than I proposed above and is
setpointer(_ref_hocvar, 'POINTER_name', nrn.Mechanism) where POINTER_name does NOT have the NMODL
SUFFIX appended. For example if the mod file has a SUFFIX foo and POINTER bar and you want the
sec(x).foo.bar POINTER to point to the hoc variable t, then one says
h.setpointer(_ref_t, 'bar', sec(x).foo)
it would be an error to use the suffix form 'bar_foo'.

For POINT_PROCESS POINTER variables the syntax is as earlier proposed:
h.setpointer(_ref_hocvar, 'POINTER_name', point_process_object)
tfoutz99
Posts: 12
Joined: Wed Dec 12, 2007 3:15 pm

Re: setpointer in Python

Post by tfoutz99 »

Thank you for this fix. It has been working for a while now. However, on long simulations, i'm getting a deallocating none fault:

Sample mod file:

Code: Select all

NEURON {
	SUFFIX test
	POINTER ex
}

PARAMETER {
}

ASSIGNED {
	ex (millivolts)
}

INITIAL {
	ex = 1
}
Sample neuron/python script

Code: Select all

from neuron import h
import sys
bob=h.Section()
bob.insert('extracellular')
bob.insert('test')
for ii in xrange(100000):
     h.setpointer(bob(0.5)._ref_e_extracellular, 'ex', bob(0.5).test);print sys.getrefcount(None)
Error:

Code: Select all

Fatal Python error: deallocating None
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: setpointer in Python

Post by hines »

Thanks and sorry. I am returning it as the value of function call and forgot to increment
its reference count. The fix has been committed to the source repository.
http://www.neuron.yale.edu/hg/neuron/nr ... a881b7692d
tfoutz99
Posts: 12
Joined: Wed Dec 12, 2007 3:15 pm

Re: setpointer in Python

Post by tfoutz99 »

Works great. Thanks!
yanou

Re: setpointer in Python

Post by yanou »

Dear All,

I am trying to make the xtra mechanism work using NEURON+Python in order to perform extracellular stimulation of a bundle of axons.
I got the xtra.mod from http://www.neuron.yale.edu/ftp/ted/neur ... nd_rec.zip.

Code: Select all

import neuron
from neuron import h
# create a section (cable)
axon = h.Section()
axon.L = 3000
axon.nseg = 100
axon.diam = 3
axon.cm = 1
axon.Ra = 200
axon.insert('hh')
axon.el_hh = -54.3
axon.gkbar_hh = 0.036
axon.insert('extracellular')
axon.insert('xtra')
In order, to be able to use xtra, I understand that I have to set the pointers im_xtra and ex_xtra. I tried to do it the two following ways, which were mentioned in this thread but I get the following error for the two
TypeError: setpointer(_ref_hocvar, 'POINTER_name', point_process or nrn.Mechanism))

Code: Select all

for sec_id in h.allsec():
   sec_id.push()
   for seg in sec_id:
       h.setpointer(seg._ref_i_membrane, 'im_xtra', seg.xtra)
       h.setpointer(seg._ref_e_extracellular, 'ex_xtra', seg.xtra)
   h.pop_section()

Code: Select all

h.setpointer(axon(0.5)._ref_i_membrane, 'im_xtra', axon(0.5).xtra)
h.setpointer(axon(0.5)._ref_e_extracellular, 'ex_xtra', axon(0.5).xtra)
Could anyone explain me why I get this error ?

I also have a second problem, if I tried to insert the mechanism directly the code of my whole model and so that I might still get after this problem, so might as well also ask for this one here. For example, I would like to drive is_xtra

Code: Select all

h.tstop = 10
h.dt = 0.0025
f = 0.1 # 100 Hz
t = np.linspace(0, h.tstop, (h.tstop+2*h.dt)/h.dt, endpoint=True)
signal = np.sin(2*math.pi*freq*t)
svec = h.Vector(self.signal)
svec.play(axon(0.5).is_xtra,t)
for this I get
AttributeError: 'nrn.Segment' object has no attribute 'is_xtra'
and if I change the last line by

Code: Select all

svec.play(axon(0.5)._ref_is_xtra,t)
I get
NameError: _ref_is_xtra was not made to point to anything at __nrnsec_0x32d86a0(0.5)
So, I am not really sure how to address the global variable is_xtra, neither if I can drive it using the function play in Python. I have seen in the documentation of xtra that it is possible to do it using hoc and I read the example in the zip folder. I also currently use the function play this way in Python but using an IClamp and it works properly.

Could you indicate me the appropriate way to use the play function in Python in order to drive the is_xtra global variable ?

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

Re: setpointer in Python

Post by ted »

xtra is a density mechanism, not a point process. To quote from Michael Hines's post above
The syntax for setting a POINTER in an NMODL density mechanism . . . is
setpointer(_ref_hocvar, 'POINTER_name', nrn.Mechanism) where POINTER_name does NOT have the NMODL
SUFFIX appended. For example if the mod file has a SUFFIX foo and POINTER bar and you want the
sec(x).foo.bar POINTER to point to the hoc variable t, then one says
h.setpointer(_ref_t, 'bar', sec(x).foo)
it would be an error to use the suffix form 'bar_foo'.
Note that your example included the suffix in the pointer name.

With regard to doing Vector play from Python, it is a mistake to try to juggle multiple chainsaws until one is well practiced at juggling a single chainsaw. Suggest you first try Vector play on a very simple model. Create a dendrite with the pas mechanism and nseg = 1, and use a pair of Vectors (a time vector and a g vector) to force g_pas at 0.5 to follow this time course:
0.001 S/cm2 for the first ms of the simulation
0.002 S/cm2 for the second ms of the simulation
and finally
0.001 S/cm2 forever after. What you learn from doing that will be very helpful in applying Vector play to more complex problems.
yanou

Re: setpointer in Python

Post by yanou »

Dear ted,

Thank you for your answer. For my first question, I should have read this thread again more carefully before asking obviously.

I have tried your little exercise because it seems there is something important to learn from it. But I think I am missing something ?

Code: Select all

import neuron
from neuron import h
h('load_file("noload.hoc")')
import numpy as np
import matplotlib.pyplot as plt

h.tstop = 1e1 # set simulation duration (ms)
h.dt = 0.0025
h.celsius = 33

axon = h.Section()
axon.L = 30
axon.nseg = 1
axon.diam = 1
axon.cm = 1
axon.Ra = 200
axon.insert('pas')
axon.e_pas = -65 #mV
axon.g_pas = 1.5e-6 # S/cm2

# record time variable
trec = h.Vector()
trec.record(h._ref_t)

g_pas_rec = h.Vector()
#g_pas_rec.record(axon(0.5)._ref_v) #  I don't know how to record g_pas variation from NEURON

t = np.linspace(0,h.tstop, round((h.tstop+h.dt)/h.dt), endpoint = True);
values = np.zeros(len(t));
values += 0.001*(t<=1);
values += 0.002*(t>1)*(t<=2);
values += 0.001 *(t>2);
svec = h.Vector(values)
svec.play(axon(0.5).g_pas,t) # equivalent to "svec.play(axon.g_pas,t)" because only one segment here.

h.run()

plt.plot(np.array(trec),svec , '-r')
plt.xlabel('Time(ms)')
plt.ylabel('g_pas (S/cm2)')
plt.figure()

plt.show()
Edit: I'm not sure but maybe the problem lies in the fact that my t vector shouldn't give all the the time value ? I mean usually in the play function I put as a second parameter h.dt like this for example: svec.play(stim._ref_amp,h.dt) where stim is h.IClamp(pos,sect). However, during this weekend I read in some threads of this forum that the second parameter was a vector, which seemed to represent the time course.

My main question would be, how can I address the global variable is of the axon section ?
-----------------------------------------
How to simulate extracellular stimulation
-----------------------------------------

Insert the extracellular and xtra mechanisms in all sections that are
subject to the extracellular field.
Compute the transfer resistance rx for every section that contains
xtra, as illustrated below.
Construct a stimulus waveform template and copy it to a Vector.
For each internal node along the axon, use this Vector to drive
is_xtra(x). The xtra mechanism uses the rx values to convert the
stimulus current waveform into the proper amplitude and sign of the
local extracellular field.

If rho, b, or c is changed, new_elec() must be invoked
to make the changes take effect.
According to the documentation of the xtra mechanism I would expect something like axon.is_xtra(0.5) but I get the following error:
AttributeError: 'nrn.Section' object has no attribute 'is_xtra'
Or something like this:

Code: Select all

svec.play(axon(0.5)._ref_is_xtra,t)
but I get:
NameError: _ref_is_xtra was not made to point to anything at __nrnsec_0x3066bd0(0.5)
which means I should look at setting the pointer for _ref_is_xtra ?
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: setpointer in Python

Post by ted »

Rather than play "guess how to do it," here's an example of Vector.play with interpolation. Be sure to read the Programmer's Reference entry on the Vector class's play method.

Code: Select all

from neuron import h, gui

dend = h.Section(name='dend') # name='dend' for nice section names in NEURON's graphs
dend.insert('pas')
# Vector.play with interpolation needs 6 points to define a rectangular pulse
# successive elements of tvec and gvec specify the x,y coordinates of the pulse waveform's "breakpoints"
# I'm sure there is a much more Pythonic way to set up these two vectors
# but the emphasis here is on results, not elegance
tvec = h.Vector(6)
gvec = h.Vector(6)
tvec.x[0] = 0
tvec.x[1] = 1
tvec.x[2] = 1
tvec.x[3] = 2
tvec.x[4] = 2
tvec.x[5] = 3
gvec.x[0] = 0.001
gvec.x[1] = 0.001
gvec.x[2] = 0.002
gvec.x[3] = 0.002
gvec.x[4] = 0.001
gvec.x[5] = 0.001

# in hoc, vector play uses a pointer (&varname) to pass values to the driven variable
# in python use _ref_varname instead
# the driven variable is dend(0.5).pas.g, so the syntax is
gvec.play(dend(0.5).pas._ref_g, tvec, 1) # nonzero third argument means play with interpolation

# no need to make things hard--use the GUI's RunControl and a Graph for quick & dirty results
h.load_file("rig.ses") # graph plots g_pas(0.5) vs. t

print "Click on the RunControl's Init & Run button"
This example involved driving the range variable g_pas. A global parameter such as the xtra mechanism's is will have the same value in all segments of all sections that have xtra, so the syntax for driving it would be
isvec.play(xtra._ref_is, tvec, 1) # nonzero third argument means play with interpolation
where isvec is the name of the vector that contains the samples of is_xtra.
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: setpointer in Python

Post by ted »

ted wrote:A global parameter such as the xtra mechanism's is will have the same value in all segments of all sections that have xtra, so the syntax for driving it would be
isvec.play(xtra._ref_is, tvec, 1) # nonzero third argument means play with interpolation
where isvec is the name of the vector that contains the samples of is_xtra.
I guess not. Wrestling with this now, will get back with something that works.
yanou

Re: setpointer in Python

Post by yanou »

Thank you ted.

I have checked the Programmer Reference regarding the play function h[url]ttps://www.neuron.yale.edu/neuron/static/new_d ... ector.play[/url] and understood that I can either use as a second parameter Dt or the vector of the time course.

What I needed was definitely the syntax to address a variable of a mechanism like "dend(0.5).pas._ref_g", I might have come across or could have induced it from my reading of The NEURON Book and different posts of the NEURON forum but didn't manage to do it. I also tried to search the forum and the documentation for such query but didn't find it. https://www.neuron.yale.edu/neuron/stat ... html#Global

Anyway, thank you.

EDIT: Just saw you gave another reply when I was posting this one. So, I will wait. I came back to this post because apparently just is can't be addressed this way in the xtra mechanism and I wanted to see if I missed something.
In [18]: axon(0.5).xtra._ref_rx
Out[18]: <hoc.HocObject at 0x35654f8>

In [19]: axon(0.5).xtra._ref_im
Out[19]: <hoc.HocObject at 0x3565d20>

In [20]: axon(0.5).xtra._ref_ex
Out[20]: <hoc.HocObject at 0x3565e88>

In [21]: axon(0.5).xtra._ref_er
Out[21]: <hoc.HocObject at 0x3580078>

In [22]: axon(0.5).xtra._ref_is
---------------------------------------------------------------------------
AttributeError Traceback (most recent call last)
<ipython-input-22-661224d9339d> in <module>()
----> 1 axon(0.5).xtra._ref_is

AttributeError: 'nrn.Mechanism' object has no attribute '_ref_is'
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: setpointer in Python

Post by ted »

The xtra mechanism's stimulus current is an NMODL GLOBAL, so its hoc name is is_xtra without the need to refer to any particular section or compartment (segment). This means it is accessible to Python using exactly the same syntax as any hoc scalar would be (see the extensive discussion of Python accessing hoc that starts with this heading
https://www.neuron.yale.edu/neuron/stat ... .HocObject
and in particular the part about
"the most common case of a hoc function holding a pointer to a variable
such as Vector.record or Vector.play. For this one needs the
_ref_varname idiom which works for any hoc variable . . . "

In other words, if isvec holds the stimulus current values and tsvec the corresponding time values,
isvec.play(h._ref_is_xtra, tsvec, 1)
will make the waveform stored jointly i nisvec and tsvec drive the xtra mechanism's is POINTER variable.

So this works

Code: Select all

# plays a waveform into xtra's is
from neuron import h, gui
dend = h.Section(name='dend') # name='dend' for nice section names in NEURON's graphs
dend.insert('pas')
dend.insert('extracellular') # this mechanism's defaults are OK
dend.insert('xtra')
dend.nseg = 5 # so we can see variation over distance
for seg in dend:
  seg.xtra.rx = seg.x*3e-5 # just for fun

# Vector.play with interpolation needs 6 points to define a rectangular pulse
# successive elements of texvec and isvec specify the x,y coordinates
# of the pulse waveform's "breakpoints"
tsvec = h.Vector(6)
isvec = h.Vector(6)
tsvec.x[0] = 0.0
tsvec.x[1] = 1.0
tsvec.x[2] = 1.0
tsvec.x[3] = 2.0
tsvec.x[4] = 2.0
tsvec.x[5] = 3.0
isvec.x[0] = 0.0
isvec.x[1] = 0.0
isvec.x[2] = 1.0
isvec.x[3] = 1.0
isvec.x[4] = 0.0
isvec.x[5] = 0.0

# xtra's ex must be coupled to e_extracellular in the same segment
# and its im must be coupled to i_membrane in the same segment
for seg in dend:
  # couple ex_xtra to e_extracellular
  # hoc syntax:  setpointer ex_xtra(x), e_extracellular(x)
  h.setpointer(dend(seg.x)._ref_e_extracellular, 'ex', dend(seg.x).xtra)
  # couple im_xtra to i_membrane
  # hoc syntax:  setpointer im_xtra(x), i_membrane(x)
  h.setpointer(dend(seg.x)._ref_i_membrane, 'im', dend(seg.x).xtra)

# now drive h.is_xtra
isvec.play(h._ref_is_xtra, tsvec, 1)

# no need to make things hard--use the GUI's RunControl and a few Graphs
# created with the NEURONMainMenu
h.load_file("rig2.ses") # graphs plot is_xtra, v, and e_extracellular vs. t
                        # at the locations of dend's internal nodes

print "Click on the RunControl's Init & Run button"
I'll email you a zip file that contains all necessary files, including rig2.ses.
yanou

Re: setpointer in Python

Post by yanou »

Thank you ted.

I have another question regarding the xtra mechanism. I don't know if it is better to create a new discussion or if I can ask it here as a follow up, which might be useful to someone looking for this mechanism in the future ? I choose the second option.

In the file calcrxc.hoc of the zip folder mentioned above there is this part of code:

Code: Select all

create sElec  // bogus section to show extracell stim/rec electrode location
proc put_sElec() {
	sElec {
		// make it 1 um long
		pt3dclear()
		pt3dadd($1-0.5, $2, $3, 1)
		pt3dadd($1+0.5, $2, $3, 1)
	}
}
put_sElec(XE, YE, ZE)

objref gElec  // will be a Shape that shows extracellular electrode location
gElec = new Shape(0)  // create it but don't map it to the screen yet
gElec.view(-245.413, -250, 520.827, 520, 629, 104, 201.6, 201.28)
objref pElec  // bogus PointProcess just to show stim location
sElec pElec = new PointProcessMark(0.5) // middle of sElec
gElec.point_mark(pElec, 2) // make it red

proc setelec() {
	xe = $1
	ye = $2
	ze = $3
	setrx(xe, ye, ze)
	put_sElec(xe, ye, ze)
}
If I understand it correctly, apart from the setrx(xe,ye,ze) call, this part seems to be to add the electrode representation in the NEURON gui.
Meaning that if I don't code the rest in Python to perform extracellular stimulation there won't be any difference ? Particularly, I don't think NEURON take into the Shape of the electrode for the different calculations of potentials and currents performed (at least not in the xtra.mod file) ?

I'm not sure of this at all, I didn't find documentation on the object PointProcessMark(), but from the name and the comments I assume it is the representation of the electrode in the GUI.

Yann
Post Reply