Numpy?

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

Moderator: hines

gartland

Numpy?

Post by gartland »

hi,
firstly, i really appreciate how easy it was for me to build my own copy of NEURON with python, install hoc as a module and get a simulation up and going from within ipython. this is essentially my first foray into ubuntu and python (with minimal hoc experience) and i've got the whole thing up and going in a couple days!

from what my coworkers using python tell me, i'm going to be wanting to do my analysis using numpy arrays (if only because the syntax and plotting will be so similar to matlab). i can't figure out a good way to export hoc vectors into a numpy array. do you have a recommendation? i did notice in a forum that you can use --enable-numpy when configuring the nrn makefile, and i've now done this, but have no idea what its good for. any tips?

thanks alot,

-Andrew
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Post by hines »

I believe that --with-numpy allows
fast copying of data between hoc Vector and Python numpy 1-d arrays. I don't have numpy on my machine and have not been maintaining that code portion so I'm not even sure it is presently working. You might inquire of Andrew Davison or Eilif Muller who worked on it about a year ago.
I'm still dithering about whether to make numpy a standard part of the distribution and will decide about it in a few months at
http://neuralensemble.org/meetings/CodeJam2.html

Til all that is straightened out, you ought to be able to plod along using the fact that if there is a hoc object which is a
Vector (or Matrix) then you can say

Code: Select all

print h.vec.size()
h.vec.x[2] = 10
print h.vec.x[2]
emuller
Posts: 15
Joined: Thu Mar 02, 2006 5:26 am
Location: Lausanne

Post by emuller »

I have just attempted a build from a recent SVN revision, and confirm that numpy support is currently broken. I will fix it and submit a patch shortly.

However, presently the following should work:

Code: Select all

import numpy
import neuron
v = neuron.Vector(10)
a = numpy.array(v)
print repr(a)
# out: array([ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.])
In the future,

Code: Select all

a = v.toarray()
Should also work, and be much faster.

I will implement v.toarray() so as to return array(v) if numpy support was not compiled into hoc.so, so one should use v.toarray() over array(v) in the near future.

More shortly.

cheers,

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

Post by hines »

I will implement v.toarray() so as to return array(v) if numpy support was not compiled into hoc.so, so one should use v.toarray() over array(v) in the near future.
That seems good, Eilif. And should we
also use v.fromarray(a) to do the inverse
(or what is the name)?
emuller
Posts: 15
Joined: Thu Mar 02, 2006 5:26 am
Location: Lausanne

Post by emuller »

I would like to point out that the class neuron.Vector already has a complete interface to play with numpy:

Code: Select all

In [1]: import neuron
NEURON -- VERSION 6.2.11 (40) 2008-02-04 (40M)
Duke, Yale, and the BlueBrain Project -- Copyright 1984-2007
See http://www.neuron.yale.edu/credits.html


In [2]: l = [1,2,3,4.5]   

In [3]: v = neuron.Vector(l)

In [4]: v
Out[4]: 1       2       3       4.5


In [5]: v[:2]
Out[5]: [1.0, 2.0]

In [6]: import numpy

In [7]: a = numpy.array([1.0,2.0,3.0])

In [8]: a[2:4]
Out[8]: array([ 3.])

In [9]: a[1:]
Out[9]: array([ 2.,  3.])

In [10]: v[:2] = a[1:]

In [11]: v
Out[11]: 2      3       3       4.5

# NEURON official release:
In [12]: v1 = neuron.Vector(list(a))
# neuralensemble.org nrnpy trunk
In [12]: v1 = neuron.Vector(a)

In [13]: v1
Out[13]: 1      2       3

In [14]: a2 = numpy.array(v1)

In [15]: a2
Out[15]: array([ 1.,  2.,  3.])

# hoc.so was built with --enable-numpy
# using the neuralensemble.org nrnpy trunk.
# It now builds, but is still broken.
In [16]: a3 = v1.toarray()
Segmentation fault
The --enable-numpy build flag is an experiment to improve performance.

The argument for why numpy support is important can be seen in the example below:

Code: Select all


In [1]: import neuron
NEURON -- VERSION 6.2.11 (40) 2008-02-04 (40M)
Duke, Yale, and the BlueBrain Project -- Copyright 1984-2007
See http://www.neuron.yale.edu/credits.html

In [2]: import numpy      

In [5]: import time

# 10 seconds at 0.01 ms resolution
In [6]: a = numpy.arange(0,10.0,0.00001)

# time is in seconds
In [7]: t1 = time.time(); b = numpy.array(a); t2 = time.time();print t2-t1
0.0187768936157

In [9]: t1 = time.time(); v = neuron.Vector(a); t2 = time.time();print t2-t1
21.2109799385

In [10]: t1 = time.time(); a = numpy.array(v); t2 = time.time();print t2-t1
43.7410931587

In [11]: numpy.alltrue(b==a)
Out[11]: True

As one sees, the conversion is way too slow.

For comparison, let's time the generation of as many random numbers as the size of a

Code: Select all


In [13]: import numpy.random

In [14]: t1 = time.time(); r = numpy.random.uniform(size=shape(a)); t2 = time.time();print t2-t1
0.112787008286

We should be able to achieve at least this performance with --enable-numpy and toarray and fromarray functions.

regards,

Eilif
Last edited by emuller on Wed Mar 19, 2008 8:05 am, edited 1 time in total.
emuller
Posts: 15
Joined: Thu Mar 02, 2006 5:26 am
Location: Lausanne

Post by emuller »

After looking into the numpy.array function, I've changed my mind regarding the creation of toarray, fromarray functions.

Code: Select all

In [17]: a = numpy.array()
In [18]: ? a.__array__
Type:           builtin_function_or_method
Base Class:     <type 'builtin_function_or_method'>
String Form:    <built-in method __array__ of numpy.ndarray object at 0xfa3100>
Namespace:      Interactive
Docstring:
    a.__array__(|dtype) -> reference if type unchanged, copy otherwise.
    
    Returns either a new reference to self if dtype is not given or a new array
    of provided data type if dtype is different from the current dtype of the
    array.


In19]: ? numpy.array

Type:           builtin_function_or_method
Base Class:     <type 'builtin_function_or_method'>
String Form:    <built-in function array>
Namespace:      Interactive
Docstring:
    array(object, dtype=None, copy=1,order=None, subok=0,ndmin=0)

    Return an array from object with the specified date-type.

    Inputs:
      object - an array, any object exposing the array interface, any
                object whose __array__ method returns an array, or any
                (nested) sequence.


This implies that if we define an __array__ method conditional on the --enable-numpy flag, we can have performance for a = numpy.array(v).

The neuron.Vector constructor should be written to automagically call a hidden .fromarray-like function in the case that the provided iterable is an array.

Therefore users can proceed as normal, and enjoy improved performance in the future with no change of model code.

cheers,

Eilif
gartland

Post by gartland »

This is all extremely helpful. I hadn't realized I could use:

Code: Select all

numpy_array=numpy.array(hoc_vec)
hoc_vec=neuron.Vector(numpy_array)
This will make my code cleaner and ready to go when you figure out a way to copy the arrays quickly! I'll just try to stick with brief models for now :)

Thanks,

-Andrew
Last edited by gartland on Thu Mar 27, 2008 4:55 pm, edited 1 time in total.
emuller
Posts: 15
Joined: Thu Mar 02, 2006 5:26 am
Location: Lausanne

Post by emuller »

The nrnpy trunk (http://neuralensemble.org/trac/nrnpy) now contains
code for optimized conversion of Hoc Vectors to numpy arrays.

Code: Select all


In [10]: a = numpy.arange(0,10.0,0.00001)

In [11]: v = neuron.Vector(a)

In [12]: t1 = time.time(); a = numpy.array(v); t2 = time.time();print t2-t1
0.0429871082306

These tests were run on the same machine as my previous post above, and so the times are directly comparable. The performance target has been surpassed.

This feature can be enabled at configure time with the flag --enable-numpy.

Please post if you have problems building.

Support for optimized conversion of numpy arrays to hoc Vectors is not yet implemented.


regards,

Eilif
gartland

Post by gartland »

I tried to build the latest version of nrnpy from your subversion repository (revision 45 i think?). I got an error when I did "make". A bunch of things appeared to "make" just fine, but this one failed and I don't have any idea why. I've had no problems building the code from the NEURON website, but haven't tried any previous revisions of the nrnpy trunk. How do I use svn to get a previous revision? I could try something earlier to see if it is unrelated to your new code (it could just be my problem right?) Let me know if I can help.

Code: Select all

Making all in e_editor
make[3]: Entering directory `/home/gartland/neuron/nrnpy/src/e_editor'
if gcc -DHAVE_CONFIG_H -I. -I. -I../.. -I../.. -I../.. -I../../src/oc -I../../src/oc -I../../src/parallel -I../../src/nrnjava -I../../src/nrncvode -I../../src/ivos -I../../src/sundials -I../../src/nrnpython -I../../src/oc    -g -O2 -MT hoc_e.o -MD -MP -MF ".deps/hoc_e.Tpo" -c -o hoc_e.o hoc_e.c; \
        then mv -f ".deps/hoc_e.Tpo" ".deps/hoc_e.Po"; else rm -f ".deps/hoc_e.Tpo"; exit 1; fi
hoc_e.c: In function ‘puts’:
hoc_e.c:2485: warning: argument ‘as’ doesn’t match built-in prototype
/bin/bash ../../libtool --tag=CC --mode=link gcc  -g -O2   -o hoc_e  hoc_e.o  -ldl -lm
mkdir .libs
gcc -g -O2 -o hoc_e hoc_e.o  -ldl -lm  
make[3]: Leaving directory `/home/gartland/neuron/nrnpy/src/e_editor'
Making all in modlunit
make[3]: Entering directory `/home/gartland/neuron/nrnpy/src/modlunit'
yacc  -d parse1.y
/bin/bash: yacc: command not found
make[3]: *** [parse1.c] Error 127
make[3]: Leaving directory `/home/gartland/neuron/nrnpy/src/modlunit'
make[2]: *** [all-recursive] Error 1
make[2]: Leaving directory `/home/gartland/neuron/nrnpy/src'
make[1]: *** [all-recursive] Error 1
make[1]: Leaving directory `/home/gartland/neuron/nrnpy'
make: *** [all] Error 2
emuller
Posts: 15
Joined: Thu Mar 02, 2006 5:26 am
Location: Lausanne

Post by emuller »

gartland wrote:

Code: Select all

Making all in e_editor
make[3]: Entering directory `/home/gartland/neuron/nrnpy/src/e_editor'
if gcc -DHAVE_CONFIG_H -I. -I. -I../.. -I../.. -I../.. -I../../src/oc -I../../src/oc -I../../src/parallel -I../../src/nrnjava -I../../src/nrncvode -I../../src/ivos -I../../src/sundials -I../../src/nrnpython -I../../src/oc    -g -O2 -MT hoc_e.o -MD -MP -MF ".deps/hoc_e.Tpo" -c -o hoc_e.o hoc_e.c; \
        then mv -f ".deps/hoc_e.Tpo" ".deps/hoc_e.Po"; else rm -f ".deps/hoc_e.Tpo"; exit 1; fi
hoc_e.c: In function ‘puts’:
hoc_e.c:2485: warning: argument ‘as’ doesn’t match built-in prototype
/bin/bash ../../libtool --tag=CC --mode=link gcc  -g -O2   -o hoc_e  hoc_e.o  -ldl -lm
mkdir .libs
gcc -g -O2 -o hoc_e hoc_e.o  -ldl -lm  
make[3]: Leaving directory `/home/gartland/neuron/nrnpy/src/e_editor'
Making all in modlunit
make[3]: Entering directory `/home/gartland/neuron/nrnpy/src/modlunit'
yacc  -d parse1.y
/bin/bash: yacc: command not found
make[3]: *** [parse1.c] Error 127
make[3]: Leaving directory `/home/gartland/neuron/nrnpy/src/modlunit'
make[2]: *** [all-recursive] Error 1
make[2]: Leaving directory `/home/gartland/neuron/nrnpy/src'
make[1]: *** [all-recursive] Error 1
make[1]: Leaving directory `/home/gartland/neuron/nrnpy'
make: *** [all] Error 2
The build error seems to be due to an incomplete build environment. Specifically, could it be you are missing yacc?

If you are using ubuntu, try installing the bison and flex packages:

Code: Select all

$ sudo apt-get install bison flex
That should fix that problem.

cheers, e.
gartland

Post by gartland »

That did the trick. Thanks!
gartland

Post by gartland »

I ran the tests that you ran in your previous post and get weird results. It seems the installation is fine (when i start ipython the version has changed...) but the function doesn't appear to work correctly.

I'm getting along okay without these fast routines for now, but I just thought I would let you know since maybe it will help you debug? I do look forward to it working though!

Code: Select all

NEURON -- VERSION 6.0.pygetsetcall.23 (35) 2007-05-23
by John W. Moore, Michael Hines, and Ted Carnevale
Duke and Yale University -- Copyright 1984-2007

Python 2.5.1 (r251:54863, Mar  7 2008, 04:10:12) 
Type "copyright", "credits" or "license" for more information.

IPython 0.8.1 -- An enhanced Interactive Python.
?       -> Introduction to IPython's features.
%magic  -> Information about IPython's 'magic' % functions.
help    -> Python's own help system.
object? -> Details about 'object'. ?object also works, ?? prints more.

  Welcome to pylab, a matplotlib-based Python environment.
  For more information, type 'help(pylab)'.

In [1]: import neuron

In [2]: import numpy

In [3]: import time

In [4]: # 10 seconds at 0.01 ms resolution

In [5]: a = numpy.arange(0,10.0,1e-5)

In [6]: # time is in seconds

In [7]: t1 = time.time(); b = numpy.array(a); t2 = time.time(); print t2-t1
0.0122640132904

In [8]: t1 = time.time(); v = neuron.Vector(a); t2 = time.time(); print t2-t1
1.4066696167e-05

In [9]: t1 = time.time(); a = numpy.array(v); t2 = time.time();print t2-t1
0.0344829559326

In [10]: numpy.alltrue(b==a)
Out[10]: False

In [11]: print b
[  0.00000000e+00   1.00000000e-05   2.00000000e-05 ...,   9.99997000e+00
   9.99998000e+00   9.99999000e+00]

In [12]: print a

/usr/lib/python2.5/site-packages/neuron/__init__.py in __getattr__(self, name)
     38 
     39     def __getattr__(self,name):
---> 40         return getattr(self.hoc_obj, name)
     41 
     42     def __len__(self):

/usr/lib/python2.5/site-packages/neuron/__init__.py in __getattr__(self, name)
     38 
     39     def __getattr__(self,name):
---> 40         return getattr(self.hoc_obj, name)
     41 
     42     def __len__(self):

#hundreds of times!

<type 'exceptions.RuntimeError'>: maximum recursion depth exceeded
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Post by hines »

Some time ago Eilif wrote that the conversion of a million length vector needs much better performance

Code: Select all

In [9]: t1 = time.time(); v = neuron.Vector(a); t2 = time.time();print t2-t1
21.2109799385 
I've been worried about making numpy a required part of the neuron distribution, e.g. in the mswin setup installer, and am doing some c implementation experiments using the Sequence protocol. Native code that knows about numpy is of course optimal but the generic sequence method is not too bad and gives performance results:

Code: Select all

a = range(0, 1000000)                         0.0476
v = neuron.Vector(a)                          0.0695
v.sum()                                   0.00307
499999500000.0
b = v.to_array()                              0.0535
len  1000000
a = numpy.arange(0, 10, .00001)               0.0132
v = neuron.Vector(a)                          0.159
v.sum()                                   0.00297
4999995.0
b = v.to_array(numpy.zeros(v.size()))         0.0699
<type 'numpy.ndarray'>
b.sum()                                   0.00875
4999995.0
This has not been committed to the subversion repository. The extension consists of two new Hoc Vector methods, vector = Vector.from_array(source)
where source should be a List or numpy 1-d array of doubles, and po = Vector.to_array() which returns a python list, or po = Vector.to_array(source) where po is the source and the source must be a python list or numpy array with the same length as the Vector.

Perhaps I should change the names to
from_python and to_python?
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Post by hines »

The subversion repository
http://www.neuron.yale.edu/cgi-bin/trac ... geset/2047
allows Vector.to_python() and Vector.from_python() methods for faster copying of python list and numpy 1-d arrays to and from hoc vectors.
gartland

Re: Numpy?

Post by gartland »

Thought I would update NEURON to stay current now that I'm getting back to the modeling again. Does this mean that the configure options "--enable-numpy" and "--with-numpy" are no longer needed or even valid?
Post Reply