Pytables/Numpy: lesson learned

Today's PyTables/NumPy lesson that I learned the hard way (i.e. through time wasted): When you use the __getitem__ method from a PyTables Table, and you pass an integer, you don't get back a record array. You get back the same thing as if you pass an integer to the __getitem__ method of a record array - namely a numpy.void instance, which is used presumably because NumPy doesn't know what to call whatever you have stashed together in one record.

Multiplying lots of matrices in NumPy

The other day I found myself needing to perform matrix multiplication in Python, using NumPy. Well, what's the big deal, you say? You do know that there exists a dot method, right?

Yes, I do know that, you smart internet person you. However, my problem was that I had a number of matrices for which I wanted to perform the same type of matrix multiplication. I had on the order of a hundred thousand five by two matrices that were to be transposed and multiplied with another hundred thousand five by two matrices.

Well, duh, you say. The dot method can handle more than two dimensions, you know. Yeah, I know that as well. However, it doesn't handle it the way I needed for this task. I wanted to end up with a hundred thousand two by two matrices. Had I used dot, I would have ended up with a hundred thousand by two by hundred thousand by two matrix.

So, I had to improvise:
>>> A.shape
(100000, 2, 5)
>>> B.shape
(100000, 2, 5)
>>> result = np.sum(np.swapaxes(A, 1, 2)[:, np.newaxis, :, :] * B[:, :, :, np.newaxis], 2)

Kind of involved, but it worked. I got the initial idea from here, but the solution given here only works for symmetrical matrices - for non-symmetrical ones you have to shift the newaxis one step to the left, or it will violate the broadcasting rules of NumPy.

Check for duplicates in Python

Today's trick: Check whether a Python container cont contains duplicates!

if len(cont) != len(set(cont)): raise myError

Neat!

NumPy structured arrays

I'm programming quite a bit in Python, and my understanding of that language is incremental. Due to the nature of my work, I also work a lot with NumPy. Today I had to solve the following problem:

• Take an input dictionary
• Create a NumPy structured array with the keys as field names, the datatypes of the values as the field datatypes, and the values themselves as the array elements.

Turns out, even though i have worked with the ndarray dtype a lot, I don't have a very good grasp of the structured array. This post will not be precise, but just a summary of what I have understood today.

First of all, I thought that a structured array would be like a 'normal' NumPy array, just that one of the dimensions had field names and data types associated with them.

But I think I am wrong in this - I think it's more a matter of a structured array being a NumPy array, where each element in the array is a structure (which makes sense once I think about it).

For instance, you can't slice a structured array according to the first interpretation:

In [1]: dtype = ''.join(('uint8,', 4*'int16,', 'int16'))
In [2]: b = np.array([(0, 1, 2, 3, 4, 5)], dtype=dtype)
In [3]: b.shape
Out[3]: (1,)
In [4]: b[0, 3]
---------------------------------------------------------------------------
IndexError Traceback (most recent call last)
----> 1 b[0, 3]
IndexError: too many indices

However, if you treat the result of a slice as a separate array, it works:

In [5]: b[0][3]
Out[5]: 3

This is very basic, I know. But it's something I learned today. And that's what this blog mainly is for. Hopefully I will learn more interesting stuff later.