Posts Tagged ‘Open Source’

GSL and C++

Tuesday, May 27th, 2008

Yesterday, in my machine learning class, I showed some example code for artificial neural networks. I am giving my students a project tomorrow where they have to do regression using neural networks, and I wanted to show them how to program these networks. The pictures we draw of networks, and the mathematics we do with them, looks less like the implementation than the other machine learning methods we cover in the course, so I thought a code demo would be in order.

The project

It's a simple feed-forward neural network I want them to implement in their project, so essentially they just need to be able to calculate target values from predictor variables and the network weights and then minimise the error function (as a function of the network weights).

target = y(x,w)
training: argminw Σn(tn - y(xn,w))2
Using the back-propagation algorithm they can get the gradient of the error function with respect to w and use that in the optimisation.Now, I could ask them to use the gradient descent algorithm or similar for the optimization, but this is not a course on numerical methods (even simple ones) and I don't want to give them the impression that numerical optimization is simpler than it really is, so I suggest to them that they use a library for the optimization.

Example implementations

For class I showed how simple neural networks could be implemented and trained in Python, using scipy to do the training. This works very fine for the demonstration, since I can interactively train the networks to various target functions and plot how well the network fits the target.

To show that the language is more or less arbitrary for the solution, I also implemented a neural network in C++, this time using GSL to minimize the error function.

In theory, there shouldn't be much difference between the Python and the C++ solutions -- and for the actual neural network implementation there isn't -- but using GSL is much more tedious than the Python approach. It really shouldn't be.

The main problem is that GSL is a C library and not a C++ library, so the features in C++ that would make the C++ solution just as simple and intuitive as the Python solution are simply not supported. One thing is operator overloading. It is missing from C, so working with vectors is a bit annoying. Especially since, in C++, you are used to the STL vectors, and switching back and forth between std::vector<> and gsl_vector * is a bit annoying. It is not the worst problem, though. Not in the neural network application, at least (it might be for more linear algebra intensive applications).

The one thing that annoyed me the most was how I should call the minimization function. In Python it was reasonably easy using high-level functions, so I could curry my error function and create a function of just the weights. In GSL, I had to go through a void pointer. In C++, I would have used a functor to deal with this, but since GSL is C it doesn't have templates, so I couldn't do that.

I wasn't happy with resource management either. GSL is clearly not designed with exceptions in mind -- after all, they do not exist in C -- and this means that you either have to wrap all GSL types in small specialized classes so they will be destroyed properly, or you need to be careful with calling the correct <type>_free() functions. The later is very hard to do exception safe.

C++ wrappers for GSL?

I googled around a bit for C++ wrappers for GSL, and it seems like there are a few around, but not really any I found particularly attractive. GSL for C++ doesn't appear particularly active (but I could be wrong). GSLwrap seems like it only wraps the linear algebra.

Are there any good GSL wrappers around, that use the high level features of C++?  Or any obvious replacement?

Just overloading operators to provide linear algebra syntax for matrix manipulation is not the way to go, unless there is some pretty heavy template magic involved, I think.  The methods in BLAS are pretty efficient because they can do a lot of operations with a single call.  Doing the same using operators would be very inefficient.  That is why the efficient C++ solutions use template expressions.  I am not sure that there is an easy way to do this through GSL, although the template expressions could end up just calling the BLAS routines, of course.

What about high level functions (functors) for optimization?

Maxima, a computer algebra system

Saturday, March 22nd, 2008

The last couple of days, I have been learning to use Maxima, a computer algebra system. I haven't really used a symbolic computer algebra system for years. I used Maple the first two years when I studied math, but after that, I never really felt that I needed it. When getting my computer science degree, that was true, but now that I am doing bioinformatics I do a lot more math than I used to.

A couple of colleagues use Mathematica, but I am too cheap for that, so I went looking for an open source alternative. I did that a few years ago, but didn't find anything I like (I even tried Maxima back then, but wasn't impressed at the time). Then, earlier this week, one of the post docs at BiRC was using wxMaxima, a GUI wrapping the Maxima system, and it looked really nice, so I decided to give it another try. This time, I'm more impressed. I'm not really sure if it is because my maths skills have improved, or just that the GUI makes it much nicer to work with, though.

Now, I don't know Mathematica and it has been ages since I used Maple, so I don't know how it compares with those. I would imagine that its feature set is limited compared to those, but I do not really need that complicated mathematics, so I doubt it is something that will be a problem for me.

I don't imagine that I will ever do any serious programming in Maxima, though. I will probably only use it to solve a few equations here and there and leave the serious data analysis to R and my numerical computation tasks to Octave or C++ with GSL.

Open data cancer

Friday, March 14th, 2008

I just saw this post on business|bytes|genes|molecyles: Open data, cancer biomarkers and a seventeen year old.

This is really cool. With open data, everyone can analyse the data and contribute to science. I have always hated that I didn't have access to data, but plenty of ideas about how to analyse it -- I am a methods development guy more than a data guy myself. These days I do have data access through a number of projects, but typically with lots of restrictions on what I can do with the data.

The more open data we have, the better...

Now I'm off to give my talk on association mapping at the department of health. I'll go searching for more info on this story later.

Clay Shirky on Love

Saturday, March 8th, 2008

A friend sent me the link for this video:

I must say, I loved it.