Mailman 3 October 2009 - NumPy-Discussion

Designing a new storage format for numpy recarrays
by Stephen Simmons 30 Oct '09

30 Oct '09

Hi, Is anyone working on alternative storage options for numpy arrays, and specifically recarrays? My main application involves processing series of large recarrays (say 1000 recarrays, each with 5M rows having 50 fields). Existing options meet some but not all of my requirements. Requirements -------------- The basic requirements are: Mandatory - fast - suitable for very large arrays (larger than can fit in memory) - compressed (to reduce disk space, read data more quickly) - seekable (can read subset of data without decompressing everything) - can append new data to an existing file - able to extract individual fields from a recarray (for when indexing or processing needs just a few fields) Nice to have - files can be split without decompressing and recompressing (e.g. distribute processing over a grid) - encryption, ideally field-level, with encryption occurring after compression - can store multiple arrays in one physical file (convenience) - portable/stardard/well documented Existing options ----------------- Over the last few years I've tried most of numpy's options for saving arrays to disk, including pickles, .npy, .npz, memmap-ed files and HDF (Pytables). None of these is perfect, although Pytables comes close: - .npy - not compressed, need to read whole array into memory - .npz - compressed but ZLIB compression is too slow - memmap - not compressed - Pytables (HDF using chunked storage for recarrays with LZO compression and shuffle filter) - can't extract individual field from a recarray - multiple dependencies (HDF, PyTables+LZO, Pyh5+LZF) - HDF is standard but LZO implementation is specific to Pytables (similarly LZF is specific to Pyh5) Are there any other options? Thoughts about a new format -------------------------------- It seems that numpy could benefit from a new storage format. My first thoughts involve: - Use chunked format - split big arrays into pages of consecutive rows, compressed separately - Get good compression ratios by shuffling data before compressing (byte 1 of all rows, then byte 2 of all rows, ...) - Get efficient access to individual fields in recarrays by compressing each recarray field's data separately (shuffling has nice side-effect of separating recarray fields' data) - Make it fast to compress and decompress by using LZO - Store pages of rows (and compressd field data within a page) using a numpy variation of IFF chunked format (e.g. used by the DjVu scanned document format version 3). For example, FORM chunk for whole file, DTYP chunk for dtype info, DIRM chunk for directory to pages holding rows, NPAG chunk for a page - The IFF structure of named chunk types allows format to be extended (other compressors than LZO, encryption, links to remote data chunks, etc) I'd appreciate any comments or suggestions before I start coding. References ----------- DjVu format - http://djvu.org/resources/ DjVu v3 format - http://djvu.org/docs/DjVu3Spec.djvu Stephen P.S. Maybe this will be too much work, and I'd be better off sticking with Pytables.....

6 7

[RFC] new function for floating point comparison
by David Cournapeau 30 Oct '09

30 Oct '09

Hi, I have added a couple of utilities for floating point comparison, to be used in unit tests mostly, and would like some comments, especially from people knowledgeable about floating point. http://github.com/cournape/numpy/tree/new_ulp_comp The main difference compared to other functions is that they are 'amplitude-independent', and use IEEE-754-specific properties. The tolerance is based on ULP, and two numbers x, y are closed depending on how many numbers are representable between x and y at the given precision. The branch contains the following new functions: * spacing(x): equivalent to the F90 intrinsic. Returns the smallest representable number needed so that spacing(x) + x > x. Spacing(1) is EPS by definition. * assert_array_almost_equal_nulp(x, y, nulp=1): assertion is defined as abs(x - y) <= nulps * spacing(max(abs(x), abs(y))). * assert_array_max_ulp(a, b, maxulp=1, dtype=None): given two numbers a and b, raise an assertion if there are more than maxulp representable numbers between a and b. They only support single and double precision - for complex number, one could arbitrarily define a distance between numbers based on nulps, say max of number of representable number for real and imag parts. Extended precision would be a bit more painful, because of the variety of implementations. I hope that they can give more robust/meaningful comparison for most of our unit tests, cheers, David

2 2

numpy loadtxt - ValueError: setting an array element with a sequence.
by TheLonelyStar 29 Oct '09

29 Oct '09

Hi, I am trying to load a tsv file using numpy.loadtxt: data = np.loadtxt('data.txt',delimiter='\t',dtype=np.float) And I get: ----------------- /usr/lib/python2.6/site-packages/numpy/lib/io.pyc in loadtxt(fname, dtype, comments, delimiter, converters, skiprows, usecols, unpack) 503 X = X.view(dtype) 504 else: --> 505 X = np.array(X, dtype) 506 507 X = np.squeeze(X) ValueError: setting an array element with a sequence. > /usr/lib/python2.6/site-packages/numpy/lib/io.py(505)loadtxt() 504 else: --> 505 X = np.array(X, dtype) 506 ---------------- I am on archlinux using 1.3.0. The file contians integers and floats sperated by tabs. Ideas? Thanks! Nathan -- View this message in context: http://www.nabble.com/numpy-loadtxt---ValueError%3A-setting-an-array-elemen… Sent from the Numpy-discussion mailing list archive at Nabble.com.

4 4

numpy loadtxt - ValueError: setting an array element with a sequence.
by Peter Schmidtke 29 Oct '09

29 Oct '09

Have you tried the numpy.fromfile function? This usually worked great for my files that had the same format than yours. ++ Peter ---------------------- PhD Student at the Molecular Modeling and Bioinformatics Group Dep. Physical Chemistry Faculty of Pharmacy University of Barcelona

2 1

Unexplained nans in matrix multiplication
by Dag Sverre Seljebotn 29 Oct '09

29 Oct '09

I'm getting (to me( very mysterious NaNs when doing matrix multiplication with certain (randomly generated) data: In [52]: a.shape, b.shape, i, j Out[52]: ((22, 1000), (1000, 22), 0, 16) In [53]: np.dot(a, b)[i,j] Out[53]: (31.322778824758661+nan*j) In [54]: np.dot(a[i,:], b[:,j]) Out[54]: (31.322778824758657+6.5017268607881213j) In [55]: np.any(np.isnan(a)), np.any(np.isnan(b)) Out[55]: (False, False) In [63]: np.max(np.abs(np.vstack((a.real, a.imag, b.real.T, b.imag.T)))) Out[63]: 4.0744710639852633 dtype is complex128. Is this a bug? Should I start looking in NumPy, ATLAS (Sage-compiled), the C compiler, the Fortran compiler...*shrug* I realize that matmul doesn't have to happen via naive vector dot products, but certainly one shouldn't hit Inf anywhere anyway? Dag Sverre

2 2

recommended way to run numpy on snow leopard
by Robin 29 Oct '09

29 Oct '09

Hi, I was wondering what the recommended way to run numpy/scipy on mac os x 10.6 is. I understood previously it was recommended to use python.org python and keep everything seperate from the system python, which worked well. But now I would like to have a 64 bit python and numpy, and there isn't one available from python.org. Also I think the python.org ones are built against the 10.4 SDK which I understand requires using gcc 4.0 - I was keen to try 4.2 to see if some of the differential performance I've seen between c and fortran goes away. So I guess the choices are system python with a virtualenv or a macports built python (or a hand built python). I'm thinking of using macports at the moment but I'm not sure how to handle preventing macports numpy from installing so I can use svn numpy. I'm not sure how the virtualenv will work with packaged installers - (ie how could I tell the wx installer to install into the virtualenv). I was wondering what others do? Cheers Robin

7 14

reading gzip compressed files using numpy.fromfile
by Peter Schmidtke 29 Oct '09

29 Oct '09

> Date: Wed, 28 Oct 2009 20:31:43 +0100 > From: Peter Schmidtke <pschmidtke(a)mmb.pcb.ub.es> > Subject: [Numpy-discussion] reading gzip compressed files using > numpy.fromfile > To: numpy-discussion(a)scipy.org > Message-ID: <fc345224bfa26132e9474287e32e083b(a)mmb.pcb.ub.es> > Content-Type: text/plain; charset="UTF-8" > > Dear Numpy Mailing List Readers, > > I have a quite simple problem, for what I did not find a solution for now. > I have a gzipped file lying around that has some numbers stored in it and I > want to read them into a numpy array as fast as possible but only a bunch > of data at a time. > So I would like to use numpys fromfile funtion. > > For now I have somehow the following code : > > > > f=gzip.open( "myfile.gz", "r" ) > xyz=npy.fromfile(f,dtype="float32",count=400) > > > So I would read 400 entries from the file, keep it open, process my data, > come back and read the next 400 entries. If I do this, numpy is complaining > that the file handle f is not a normal file handle : > OError: first argument must be an open file > > but in fact it is a zlib file handle. But gzip gives access to the normal > filehandle through f.fileobj. > > So I tried xyz=npy.fromfile(f.fileobj,dtype="float32",count=400) > > But there I get just meaningless values (not the actual data) and when I > specify the sep=" " argument for npy.fromfile I get just .1 and nothing > else. > > Can you tell me why and how to fix this problem? I know that I could read > everything to memory, but these files are rather big, so I simply have to > avoid this. > > Thanks in advance. > > > -- > > Peter Schmidtke > > ---------------------- > PhD Student at the Molecular Modeling and Bioinformatics Group > Dep. Physical Chemistry > Faculty of Pharmacy > University of Barcelona > > > > ------------------------------ > > Message: 2 > Date: Wed, 28 Oct 2009 14:33:11 -0500 > From: Robert Kern <robert.kern(a)gmail.com> > Subject: Re: [Numpy-discussion] reading gzip compressed files using > numpy.fromfile > To: Discussion of Numerical Python <numpy-discussion(a)scipy.org> > Message-ID: > <3d375d730910281233r5cadd0fcubea14676a3a978f1(a)mail.gmail.com> > Content-Type: text/plain; charset=UTF-8 > > On Wed, Oct 28, 2009 at 14:31, Peter Schmidtke <pschmidtke(a)mmb.pcb.ub.es> > wrote: >> Dear Numpy Mailing List Readers, >> >> I have a quite simple problem, for what I did not find a solution for >> now. >> I have a gzipped file lying around that has some numbers stored in it and >> I >> want to read them into a numpy array as fast as possible but only a bunch >> of data at a time. >> So I would like to use numpys fromfile funtion. >> >> For now I have somehow the following code : >> >> >> >> ? ? ? ?f=gzip.open( "myfile.gz", "r" ) >> xyz=npy.fromfile(f,dtype="float32",count=400) >> >> >> So I would read 400 entries from the file, keep it open, process my data, >> come back and read the next 400 entries. If I do this, numpy is >> complaining >> that the file handle f is not a normal file handle : >> OError: first argument must be an open file >> >> but in fact it is a zlib file handle. But gzip gives access to the normal >> filehandle through f.fileobj. > > np.fromfile() requires a true file object, not just a file-like > object. np.fromfile() works by grabbing the FILE* pointer underneath > and using C system calls to read the data, not by calling the .read() > method. > >> So I tried ?xyz=npy.fromfile(f.fileobj,dtype="float32",count=400) >> >> But there I get just meaningless values (not the actual data) and when I >> specify the sep=" " argument for npy.fromfile I get just .1 and nothing >> else. > > This is reading the compressed data, not the data that you want. > >> Can you tell me why and how to fix this problem? I know that I could read >> everything to memory, but these files are rather big, so I simply have to >> avoid this. > > Read in reasonably-sized chunks of bytes at a time, and use > np.fromstring() to create arrays from them. > > -- > Robert Kern > > "I have come to believe that the whole world is an enigma, a harmless > enigma that is made terrible by our own mad attempt to interpret it as > though it had an underlying truth." > -- Umberto Eco > > > ------------------------------ > > Message: 3 > Date: Wed, 28 Oct 2009 13:26:41 -0700 > From: Christopher Barker <Chris.Barker(a)noaa.gov> > Subject: Re: [Numpy-discussion] reading gzip compressed files using > numpy.fromfile > To: Discussion of Numerical Python <numpy-discussion(a)scipy.org> > Message-ID: <4AE8A901.3060403(a)noaa.gov> > Content-Type: text/plain; charset=UTF-8; format=flowed > > Robert Kern wrote: >>> f=gzip.open( "myfile.gz", "r" ) >>> xyz=npy.fromfile(f,dtype="float32",count=400) > >> Read in reasonably-sized chunks of bytes at a time, and use >> np.fromstring() to create arrays from them. > > Something like: > > count = 400 > xyz = np.fromstring(f.read(count*4), dtype=np.float32) > > should work (untested...) > > -Chris > > > > > -- > Christopher Barker, Ph.D. > Oceanographer > > Emergency Response Division > NOAA/NOS/OR&R (206) 526-6959 voice > 7600 Sand Point Way NE (206) 526-6329 fax > Seattle, WA 98115 (206) 526-6317 main reception > > Chris.Barker(a)noaa.gov > > Thanks Robert and Chris...indeed I managed to read it quite fast this way. ++ Peter Schmidtke ---------------------- PhD Student at the Molecular Modeling and Bioinformatics Group Dep. Physical Chemistry Faculty of Pharmacy University of Barcelona

1 0

Segfault when using scipy.special.hermite?
by Ole Streicher 28 Oct '09

28 Oct '09

Hi, Is there something wrong with scipy.special.hermite? The following code produces glibc errors: ------------8<----------------------- import scipy.special h = [] for i in xrange(15): print i h.append(scipy.special.hermite(i+1)) ------------8<----------------------- results in ... 12 *** glibc detected *** python: free(): invalid next size (fast): 0x00000000007e2290 *** OS: OpenSUSE 11.1 (x86_64) Python 2.6.0 Scipy: 0.7.0 When using ipython 0.8.4 on the same machine, the error does not occur. What may be the problem here? Regards Ole

3 2

reading gzip compressed files using numpy.fromfile
by Peter Schmidtke 28 Oct '09

28 Oct '09

Dear Numpy Mailing List Readers, I have a quite simple problem, for what I did not find a solution for now. I have a gzipped file lying around that has some numbers stored in it and I want to read them into a numpy array as fast as possible but only a bunch of data at a time. So I would like to use numpys fromfile funtion. For now I have somehow the following code : f=gzip.open( "myfile.gz", "r" ) xyz=npy.fromfile(f,dtype="float32",count=400) So I would read 400 entries from the file, keep it open, process my data, come back and read the next 400 entries. If I do this, numpy is complaining that the file handle f is not a normal file handle : OError: first argument must be an open file but in fact it is a zlib file handle. But gzip gives access to the normal filehandle through f.fileobj. So I tried xyz=npy.fromfile(f.fileobj,dtype="float32",count=400) But there I get just meaningless values (not the actual data) and when I specify the sep=" " argument for npy.fromfile I get just .1 and nothing else. Can you tell me why and how to fix this problem? I know that I could read everything to memory, but these files are rather big, so I simply have to avoid this. Thanks in advance. -- Peter Schmidtke ---------------------- PhD Student at the Molecular Modeling and Bioinformatics Group Dep. Physical Chemistry Faculty of Pharmacy University of Barcelona

3 2

Using matplotlib's prctile on masked arrays
by Gökhan Sever 28 Oct '09

28 Oct '09

Hello, Consider this sample two columns of data: 999999.9999 999999.9999 999999.9999 999999.9999 999999.9999 999999.9999 999999.9999 1693.9069 999999.9999 1676.1059 999999.9999 1621.5875 651.8040 1542.1373 691.0138 1650.4214 678.5558 1710.7311 621.5777 999999.9999 644.8341 999999.9999 696.2080 999999.9999 Putting into this data into a file say "sample.data" and loading with: a,b = np.loadtxt('sample.data', dtype="float").T I[16]: a O[16]: array([ 1.00000000e+06, 1.00000000e+06, 1.00000000e+06, 1.00000000e+06, 1.00000000e+06, 1.00000000e+06, 6.51804000e+02, 6.91013800e+02, 6.78555800e+02, 6.21577700e+02, 6.44834100e+02, 6.96208000e+02]) I[17]: b O[17]: array([ 999999.9999, 999999.9999, 999999.9999, 1693.9069, 1676.1059, 1621.5875, 1542.1373, 1650.4214, 1710.7311, 999999.9999, 999999.9999, 999999.9999]) ### interestingly, the second column is loaded as it is but a values reformed a little. Why this could be happening? Any idea? Anyways, back to masked arrays: I[24]: am = ma.masked_values(a, value=999999.9999) I[25]: am O[25]: masked_array(data = [-- -- -- -- -- -- 651.804 691.0138 678.5558 621.5777 644.8341 696.208], mask = [ True True True True True True False False False False False False], fill_value = 999999.9999) I[30]: bm = ma.masked_values(b, value=999999.9999) I[31]: am O[31]: masked_array(data = [-- -- -- -- -- -- 651.804 691.0138 678.5558 621.5777 644.8341 696.208], mask = [ True True True True True True False False False False False False], fill_value = 999999.9999) So far so good. A few basic checks: I[33]: am/bm O[33]: masked_array(data = [-- -- -- -- -- -- 0.422662755126 0.418689311712 0.39664667346 -- -- --], mask = [ True True True True True True False False False True True True], fill_value = 999999.9999) I[34]: mean(am/bm) O[34]: 0.41266624676580849 Unfortunately, matplotlib.mlab's prctile cannot handle this division: I[54]: prctile(am/bm, p=[5,25,50,75,95]) O[54]: array([ 3.96646673e-01, 6.21577700e+02, 1.00000000e+06, 1.00000000e+06, 1.00000000e+06]) This also results with wrong looking box-and-whisker plots. Testing further with scipy.stats functions yields expected correct results: I[55]: stats.scoreatpercentile(am/bm, per=5) O[55]: 0.40877012449846228 I[49]: stats.scoreatpercentile(am/bm, per=25) O[49]: masked_array(data = --, mask = True, fill_value = 1e+20) I[56]: stats.scoreatpercentile(am/bm, per=95) O[56]: masked_array(data = --, mask = True, fill_value = 1e+20) Any confirmation? -- Gökhan

3 5

NumPy-Discussion October 2009