Mailman 3 October 2010 - NumPy-Discussion

SVD problem - matrices are not aligned
by Daniel Wagner Oct. 24, 2010

Oct. 24, 2010

Hi, I'm a new subscriber of this list. I hope to directly start with a question is ok... My question or problem: I've a matrix A which is calculated from the data b. The shapes of these matrices are: >>>A.shape (954, 9) >>>b.shape (954,) I calculate the SVD of A: >>> U, w, V = numpy.linalg.svd(A, full_matrices="True") >>>U.shape (954, 954) >>>W.diag(w) >>>W.shape (9,9) >>>V.shape (9,9) If I'm doing the check of the SVD results using: >>>numpy.allclose(A, numpy.dot(U, numpy.dot(W, V))) I get this error: "ValueError: matrices are not aligned" I'just don't get where the error does come from. There seems to be a general problem with my data or how I'm using the svd function. I hope anybody has an idea what I'm doing wrong or where my problem comes from. I really appreciate any suggestions! Many thanks! Greetings, Daniel Wagner

3 9

Re: [Numpy-discussion] numpy.polynomial.Polynomial
by josef.pktd＠gmail.com Oct. 22, 2010

Oct. 22, 2010

trying again, my reply got bounced back by postfix program if the mail delivery service On Fri, Oct 22, 2010 at 2:45 PM, <josef.pktd(a)gmail.com> wrote: > On Fri, Oct 22, 2010 at 2:14 PM, Charles R Harris > <charlesr.harris(a)gmail.com> wrote: >> >> >> On Fri, Oct 22, 2010 at 11:47 AM, <josef.pktd(a)gmail.com> wrote: >>> >>> On Fri, Oct 22, 2010 at 12:26 PM, Charles R Harris >>> <charlesr.harris(a)gmail.com> wrote: >>> > >>> > >>> > On Fri, Oct 22, 2010 at 9:51 AM, <josef.pktd(a)gmail.com> wrote: >>> >> >>> >> I'm subclassing numpy.polynomial.Polynomial. So far it works well. >>> >> >>> >> One question on inplace changes >>> >> >>> >> Is it safe to change coef directly without creating a new instance? >>> >> I'm not trying to change anything else in the polynomial, just for >>> >> example pad, truncate or invert the coef inplace, e.g >>> >> >>> >> def pad(self, maxlag): >>> >> self.coef = np.r_[self.coef, np.zeros(maxlag - len(self.coef))] >>> >> >>> >> Currently, I have rewritten this to return a new instance. >>> >> >>> > >>> > You can (currently) modify the coef and it should work, but I think it >>> > best >>> > to regard the Polynomial class as immutable. I'm even contemplating >>> > making >>> > the coef attribute read only just to avoid such things. Another tip is >>> > to >>> > use // instead of / for division, polynomials are rather like integers >>> > that >>> > way and don't have a true divide so plain old / will fail for python 3.x >>> > >>> > Note that most operations will trim trailing zeros off the result. >>> > >>> > In [6]: P((1,1,1,0,0,0)) >>> > Out[6]: Polynomial([ 1., 1., 1., 0., 0., 0.], [-1., 1.]) >>> > >>> > In [7]: P((1,1,1,0,0,0)) + 1 >>> > Out[7]: Polynomial([ 2., 1., 1.], [-1., 1.]) >>> > >>> > The reason the constructor doesn't was because trailing zeros can be of >>> > interest in least squares fits. Is there a particular use case for which >>> > trailing zeros are important for you? The polynomial modules aren't >>> > finished >>> > products yet, I can still add some functionality if you think it useful. >>> >>> I need "long" division, example was A(L)/B(L) for lag-polynomials as I >>> showed before. >>> >> >> OK, I kinda thought that was what you wanted. It would be a version of >> "true" division, the missing pieces are how to extend that to other basis, >> there are several possibilities... But I suppose they could just be marked >> not implemented for the time being. There also needs to be a way to specify >> "precision" and the location of the "decimal" point. > As long as subclassing works and it seems to work well so far, adding a few topic specific methods is quite easy. > >> >>> >>> My current version (unfinished since I got distracted by >>> stats.distribution problems): >>> >>> from numpy import polynomial as npp >>> >>> >>> class LagPolynomial(npp.Polynomial): >>> >>> #def __init__(self, maxlag): >>> >>> def pad(self, maxlag): >>> return LagPolynomial(np.r_[self.coef, >>> np.zeros(maxlag-len(self.coef))]) >>> >>> def padflip(self, maxlag): >>> return LagPolynomial(np.r_[self.coef, >>> np.zeros(maxlag-len(self.coef))][::-1]) >>> >>> def flip(self): >>> '''reverse polynomial coefficients >>> ''' >>> return LagPolynomial(self.coef[::-1]) >>> >>> def div(self, other, maxlag=None): >>> '''padded division, pads numerator with zeros to maxlag >>> ''' >>> if maxlag is None: >>> maxlag = max(len(self.coef), len(other.coef)) + 1 >>> return (self.padflip(maxlag) / other.flip()).flip() >>> >>> def filter(self, arr): >>> return (self * arr) #trim to end >>> >>> >>> another method I haven't copied over yet is the adjusted fromroots >>> (normalized lag-polynomial from roots) >>> >>> Essentially, I want to do get the AR and ARMA processes in several >>> different ways because I don't trust (my interpretation) of any single >>> implementation and eventually to see which one is fastest. >>> >> >> >> I could also implement "polyz" polynomials that would use negative powers of >> z. The Chebyshev polynomials are currently implemented with symmetric >> z-series using both positive and negative powers, but I may change that. My background for this is pretty much causal filters in time series analysis. I still have only vague ideas about some of the signaling and polynomial stuff discussed in the previous thread. But I take whatever I can get, and can figure out how to use it. The polynomial class (and my wrappers around scipy.signal and fft) is nice because it allows almost literal translation of textbook formulas. If I have enough time, spectral densities are one of the next on the schedule. Thanks, I will keep treating the Polynomials as immutable. Josef >> >> Another possibility is some sort of factory function that emits polynomial >> classes with certain additional properties, I'm thinking of something like >> that for Jacobi polynomials. >> >> Chuck >> >> _______________________________________________ >> NumPy-Discussion mailing list >> NumPy-Discussion(a)scipy.org >> http://mail.scipy.org/mailman/listinfo/numpy-discussion >> >> >

1 1

numpy.polynomial.Polynomial
by josef.pktd＠gmail.com Oct. 22, 2010

Oct. 22, 2010

I'm subclassing numpy.polynomial.Polynomial. So far it works well. One question on inplace changes Is it safe to change coef directly without creating a new instance? I'm not trying to change anything else in the polynomial, just for example pad, truncate or invert the coef inplace, e.g def pad(self, maxlag): self.coef = np.r_[self.coef, np.zeros(maxlag - len(self.coef))] Currently, I have rewritten this to return a new instance. Thanks, Josef

2 3

dtype subarray comparison
by Mark Wiebe Oct. 21, 2010

Oct. 21, 2010

It turns out that when comparing dtypes, the subarray is currently ignored. This means you can get things like this: >>> import numpy as np >>> np.dtype(('f4',(6))) == np.dtype(('f4',(2,3))) True >>> np.dtype(([('a','i4')],2)) == np.dtype(([('b','u2'),('c','i2')],2)) True which pretty clearly should be False in both cases. I've implemented a patch to fix this, for which a pull request is here: http://github.com/numpy/numpy/pull/6 The main points included in this patch are: * If there's a subarray shape, it has to match exactly. * Currently, sometimes the subarray shape was an integer, other times a tuple. The output formatting code checks for this and prints a tuple in both cases. I changed the code to turn the integer into a tuple on construction instead. This didn't cause any tests to fail, so I think it's a safe change. * The current code converts (type, 1) and (type, tuple()) to just type, so this means (type, 1) != (type, (1,)) but (type, 2) == (type, (2,)) I put in some tests of these edge cases. http://github.com/numpy/numpy/blob/master/numpy/core/src/multiarray/descrip… Cheers, Mark

3 10

test_creation Errors with np 1.5.1rc1 built on osx 10.5 py27 32bit
by Vincent Davis Oct. 21, 2010

Oct. 21, 2010

Build info numpy-1.5.1rc1-py2.7-python.org-macosx10.5 python-2.7-macosx10.3 os x 10.5 installed on osx 10.6 python-2.7-macosx10.3 Been working on building the dmg, you can try it/download at http://vincentdavis.info/Shared/numpy-1.5.1rc1-py2.6-python.org-macosx10.3.… Saw a few tickets in trac but have not look into this really. I am thought maybe someone might already know what this is. Running unit tests for numpy NumPy version 1.5.1rc1 NumPy is installed in /Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy Python version 2.7 (r27:82508, Jul 3 2010, 20:17:05) [GCC 4.0.1 (Apple Inc. build 5493)] nose version 0.11.2 ...EEEEEEEEEEEE......................................................................................................................................................................................................................................................................................................................................................................................................................................................................................K........................................................................................................................................................................................................................................K......................Warning: overflow encountered in ldexp .......................................................................K......................K........................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................S............................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................. ====================================================================== ERROR: test_creation (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 10, in test_creation dt1 = np.dtype('M8[750%s]'%unit) TypeError: data type not understood ====================================================================== ERROR: test_creation_overflow (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 62, in test_creation_overflow timesteps = np.array([date], dtype='datetime64[s]')[0].astype(np.int64) TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_as (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 58, in test_divisor_conversion_as self.assertRaises(ValueError, lambda : np.dtype('M8[as/10]')) File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/unittest/case.py", line 456, in assertRaises callableObj(*args, **kwargs) File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 58, in <lambda> self.assertRaises(ValueError, lambda : np.dtype('M8[as/10]')) TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_bday (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 32, in test_divisor_conversion_bday assert np.dtype('M8[B/12]') == np.dtype('M8[2h]') TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_day (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 37, in test_divisor_conversion_day assert np.dtype('M8[D/12]') == np.dtype('M8[2h]') TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_fs (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 54, in test_divisor_conversion_fs assert np.dtype('M8[fs/100]') == np.dtype('M8[10as]') TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_hour (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 42, in test_divisor_conversion_hour assert np.dtype('m8[h/30]') == np.dtype('m8[2m]') TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_minute (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 46, in test_divisor_conversion_minute assert np.dtype('m8[m/30]') == np.dtype('m8[2s]') TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_month (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 21, in test_divisor_conversion_month assert np.dtype('M8[M/2]') == np.dtype('M8[2W]') TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_second (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 50, in test_divisor_conversion_second assert np.dtype('m8[s/100]') == np.dtype('m8[10ms]') TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_week (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 26, in test_divisor_conversion_week assert np.dtype('m8[W/5]') == np.dtype('m8[B]') TypeError: data type not understood ====================================================================== ERROR: test_divisor_conversion_year (test_datetime.TestDateTime) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/tests/test_datetime.py", line 16, in test_divisor_conversion_year assert np.dtype('M8[Y/4]') == np.dtype('M8[3M]') TypeError: data type not understood ---------------------------------------------------------------------- Ran 3018 tests in 10.218s FAILED (KNOWNFAIL=4, SKIP=1, errors=12) <nose.result.TextTestResult run=3018 errors=12 failures=0> -- Thanks Vincent Davis 720-301-3003

1 1

function name as parameter
by Thomas Kirk Gamble Oct. 21, 2010

Oct. 21, 2010

I'm trying to write an implementation of the amoeba function from numerical recipes and need to be able to pass a function name and parameter list to be called from within the amoeba function. Simply passing the name as a string doesn't work since python doesn't know it is a function and throws a typeerror. Is there something similar to IDL's 'call_function' routine in python/numpy or a pythonic/numpy means of passing function names? -- Thomas K. Gamble Research Technologist, System/Network Administrator Chemical Diagnostics and Engineering (C-CDE) Los Alamos National Laboratory MS-E543,p:505-665-4323 f:505-665-4267 There cannot be a crisis next week. My schedule is already full. Henry Kissinger

8 15

Count occurrences in a nested ndarray
by Dewald Pieterse Oct. 20, 2010

Oct. 20, 2010

Hi All I have a nested array created using: edges = scipy.misc.pilutil.imfilter(prent,'find_edges') edges array([[[ 0, 255, 0], [ 0, 255, 0], [ 0, 255, 0], ..., [ 0, 255, 0], [ 0, 255, 0], [ 0, 255, 0]], [[ 0, 255, 0], [ 0, 0, 0], [ 0, 0, 0], ..., [ 0, 0, 0], [ 0, 0, 0], [ 0, 255, 0]], [[ 0, 255, 0], [ 0, 0, 0], [ 0, 0, 0], ..., [ 0, 0, 0], [ 0, 0, 0], [ 0, 255, 0]], ..., [[ 0, 255, 0], [ 0, 0, 0], [ 0, 0, 0], ..., [ 0, 0, 0], [ 0, 0, 0], [ 0, 255, 0]], [[ 0, 255, 0], [ 0, 0, 0], [ 0, 0, 0], ..., [ 0, 0, 0], [ 0, 0, 0], [ 0, 255, 0]], [[ 0, 255, 0], [ 0, 255, 0], [ 0, 255, 0], ..., [ 0, 255, 0], [ 0, 255, 0], [ 0, 255, 0]]], dtype=uint8) I want to count the number of occurrences of certain unique elements in the array, I know what the elements are that I want to count [0,255,0], [255,0,0] and [0,0,255]. I want to count the number of pixels of a particular color in a picture to determine the edge length and calculate areas and such. array does not have a array.count() method and trying to count the instances using a nested for loop like: for xiter in range(xindex): for yiter in range(yindex): if edges[xiter,yiter,:] == [255,0,0]: groenpixelarea = groenpixelarea + 1 if edges[xiter,yiter,:] == [0,255,0]: rooipixelarea = rooipixelarea + 1 if edges[xiter,yiter,:] == [0,0,255]: bloupixelarea = bloupixelarea + 1 results in: 16 for xiter in range(xindex): 17 for yiter in range(yindex): ---> 18 if edges[xiter,yiter,:].any() == [255,0,0]: 19 groenpixelarea = groenpixelarea + 1 20 if edges[xiter,yiter,:] == [0,255,0]: ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all() WARNING: Failure executing file: <analiseerverwerkteprent.py> What am I doing wrong? Thanx -- Dewald Pieterse

3 3

numpy.test() segfaults '2.0.0.dev8716'
by Nils Wagner Oct. 20, 2010

Oct. 20, 2010

test_duplicate_field_names_assign (test_regression.TestRegression) ... ........FF............ Program received signal SIGSEGV, Segmentation fault. 0x00007ffff7b3077a in ?? () from /usr/lib64/libpython2.6.so.1.0 (gdb) bt #0 0x00007ffff7b3077a in ?? () from /usr/lib64/libpython2.6.so.1.0 #1 0x00007ffff7b2ea3f in _PyArg_ParseTuple_SizeT () from /usr/lib64/libpython2.6.so.1.0 #2 0x00007ffff62ce02e in VOID_getitem (ip=0x1583150 "", ap=0x175f160) at numpy/core/src/multiarray/arraytypes.c.src:565 #3 0x00007ffff62e64eb in array_toscalar (self=0x175f160, args=<value optimized out>) at numpy/core/src/multiarray/methods.c:554 #4 0x00007ffff7b2165c in PyEval_EvalFrameEx () from /usr/lib64/libpython2.6.so.1.0 #5 0x00007ffff7b26251 in PyEval_EvalCodeEx () from /usr/lib64/libpython2.6.so.1.0 #6 0x00007ffff7b2103d in PyEval_EvalFrameEx () from /usr/lib64/libpython2.6.so.1.0 #7 0x00007ffff7b26431 in PyEval_EvalCodeEx () from /usr/lib64/libpython2.6.so.1.0 #8 0x00007ffff7af3db2 in ?? () from /usr/lib64/libpython2.6.so.1.0 #9 0x00007ffff7ae7cb2 in PyObject_Call () from /usr/lib64/libpython2.6.so.1.0 #10 0x00007ffff7b1fd26 in PyEval_CallObjectWithKeywords () from /usr/lib64/libpython2.6.so.1.0 #11 0x00007ffff62c2184 in array_str (self=<value optimized out>) at numpy/core/src/multiarray/arrayobject.c:393 #12 0x00007ffff7b002e8 in _PyObject_Str () from /usr/lib64/libpython2.6.so.1.0 #13 0x00007ffff7b003ea in PyObject_Str () from /usr/lib64/libpython2.6.so.1.0 ---Type <return> to continue, or q <return> to quit---q Quit

2 1

datarray lightning talk at NYC Python group
by M Trumpis Oct. 20, 2010

Oct. 20, 2010

Just a heads up that I'm going to do a quick pitch for datarray this Wed. at the NYC Python Meetup thing. If any collaborator is attending, please introduce yourself! Mike

3 3

ANN: NUmexpr 1.4.1 released
by Francesc Alted Oct. 20, 2010

Oct. 20, 2010

========================== Announcing Numexpr 1.4.1 ========================== Numexpr is a fast numerical expression evaluator for NumPy. With it, expressions that operate on arrays (like "3*a+4*b") are accelerated and use less memory than doing the same calculation in Python. What's new ========== This is a maintenance release. On it, several improvements have been done in order to prevent deadlocks in new threaded code (fixes #33). Also the GIL is released now during computations, which should be interesting for embedding numexpr in threaded Python apps. In case you want to know more in detail what has changed in this version, see: http://code.google.com/p/numexpr/wiki/ReleaseNotes or have a look at RELEASE_NOTES.txt in the tarball. Where I can find Numexpr? ========================= The project is hosted at Google code in: http://code.google.com/p/numexpr/ And you can get the packages from PyPI as well: http://pypi.python.org/pypi Share your experience ===================== Let us know of any bugs, suggestions, gripes, kudos, etc. you may have. Enjoy! -- Francesc Alted

1 0