A Programming Language

Since I'm very interested in numpy and vectorization, I became curious about when how this approach appeared. I find this notion very convenient and very producing, while it leaves the space for optimizations inside.

Obviously, a predecessor of numpy is MATLAB, which appeared a long time ago - in the late 70's, while the first release happened more than 30 years ago - in 1984! Many things were inherited from this language (not the multiCPU/GPU backend, unfortunately, but I hope that things will change oce).

But more interesting moment that vector-based syntax was not fortran, as I thought earlier. The latter got it's operations over vectors only in Fortran'90 (which was much later). 

The real predecessor of MATLAB was APL (A Programming Language, not named after Apple company), this language was quite complicated and short,  and mostly this was looking not like a code, but as numerous formulas. Apart from this,  APL used many specific symbols and combined operators. For instance,  +\ seq is sum of elements in sequence.

Wikipedia claims that this returns list of prime numbers in $1, 2, ..R$. 

(~R∊R∘.×R)/R←1↓ιR

Now you understand why this language isn't popular today :)

From that moment we learnt that many things like map, where, sum can be written using words, and this will be much more clear and reliable.

Thing I also find interesting is APL interpretation process was much more complicated then one of numpy or matlab, and in particular supported lazy evaluation (thus, construction of expressions was available). This also means that different optimization were possible, like usage of multiple different vector operations at once (like a + b*c), if processor supports this. In russian this is called векторное зацепление, but I haven't found appropriate translation in wiki.

Numpy tricks (again)

A small post about different tricks I use in numpy. Since frequently I work with statistical data and different statistical criterions, frequently I have to compute the order statistics (alternatively, one can ask about which position will get every number in array after array was sorted)
There is rankdata function inside scipy which solves the named problem:

result = rankdata(array)

However, one can easily implement this is numpy (isn't beautiful?):

result = numpy.argsort(numpy.argsort(array))

which is (almost) equivalent code.
See also my previous post on numpy tricks
PS. Once upon a time I will understand why all blogging platforms suggest heavy and ugly layout templates. Even adding background image turns into non-trivial task.

Code and the importance of vectorization

That awkward moment, when the code written in matlab is easier to read and understand then tons os explanations:
Salakhudinov's code of RBM

This code IMHO is a good argument when you need to explain someone that he really needs to learn at least one language with vectorization, no matter whether it is R, matlab, or numpy or theano in python.

Numpy exercises

When one starts writing in python, the typical reaction is disappointment about how slow it is compared to any compilable language. After a while, you learn numpy and find out it's actually not so bad.

Having spent a month with numpy, I found out that many things can be written in it.
Having spent a year with it, I found out that almost any algorithm may be vectorized, though it's sometimes non-trivial.

I'm still quite disappointed about the most answers at stackoverflow, where people prefer plain python for any nontrivial thing.

For instance, you need to compute statistics of values in array. Well, you can sort array and keep track of initial position. Alternatively, you can do it in numpy-one-liner:

order_statistics = numpy.argsort(numpy.argsort(initial_array))

Don't believe? Check this!
Want to compute mean value over the group of events? Another one-liner:

means = numpy.bincount(group_indices, weights=values) / numpy.bincount(group_indices)

Writing oblivious decision tree in numpy is very simple and can be done really fast.

As a non-trivial problem: will you be able to write application of usual decision tree in pure numpy? For simplicity, you can first consider only trees with equal depth of all leaves.

Optimization of vector operations

Recently worked over optimization of some (secret) classifier. The problem was mostly not in training, but in applying of trained classifier — this code was originally written in C++ and then translated to cython (which decreased the speed by factor of 2).

This was quite easy rewrite the code using numpy and vectorized approach (initally predictions were built event-by event, after rewriting the classifier was applied tree-by-tree). However this gave only speed comparable with original C++ code (so, twixe faster cython).

What really fastened the code is switching from int8 operations to int64 (the latest are natively supported in all modern processors). So 8 operations in int8 were grouped into one 64-bit. This is done pretty simple by creating views:

In[1]: import numpy
In[2]: x = numpy.random.randint(0, 100, size=64000).astype(’int8′)
In[3]: y = numpy.random.randint(0, 100, size=64000).astype(’int8′)
In[4]: %timeit x & y 
10000 loops, best of 3: 60.6 µs per loop
In[5]: %timeit x.view(’int64′) & y.view(’int64′) 
100000 loops, best of 3: 12 µs per loop
# Checked that output is correct
In[6]: numpy.all( (x & y) == (x.view(’int64′) & y.view(’int64’)).view(’int8′) )
Out[6]: True

 

In this simple example we see 5x speed up. Views of course do not copy the data, which is very essential for the speed. This trick can be applied with: summation / subtration / binary or / binary and, but you need that the size of original array was divisible by 8.

By the way, there is awesome collection of twiddling bits, which was my starting point in bit optimizations.