numpy/scipy

The more I use numpy and scipy the more I like them. Matplotlib has some quirks, but all in all its easy to quickly write something powerful that just works.

Update (09 July 2009): That said, I am surprised at how simplistic the fitting routines (aka optimization, minimization, etc) are in scipy, and the lack of good generalized fitting packages in python. Its all very well if you want to use the routines provided in scipy.optimize and uncertainties in your data aren't a big issue but that's not at all the realm of data I experience daily.

It may be that scipy.optimize's routines can handle data uncertainties properly through some more complicated interface, but all the examples I've seen don't seem to suggest that as they all totally ignore data errors. Nor do they mention even simple fitting concepts like chi squared (I'd survive with chi squared, I don't need log likelihood all the time).

Update (10 July 2009): OK, you can handle data uncertainties using the scipy.optimize routines, you just make sure the function you want to minimize calculates chi squared or log likelihood, and you effectively have chi^2 fitting or Cash statistic fitting. I haven't yet worked out how to calculate the uncertainties on the fit parameters using scipy.optimize -- the most naive way is to use the covariance matrix (not robust, I know).

Space Telescope's stsci_python provides pytools.nmpfit (a numpy version of python mpfit, itself a version of IDL mpfit). This does Levenberg-Marquardt fitting and does provide estimates of the best-fit parameter uncertainties (using the covariance matrix). After giving it a spin I find it works just as well as scipy.optimize, although at the cost of having to install the entire stsci_python package. More interestingly (n)mpfit has a well developed interface for fixing (and unfixing) fit parameters as well as handling limits on the fit parameters. I ended up using this for the HST ACS data analysis routines I'm writing. Still, none of these approach the sophistication of xspec or sherpa. If it weren't for the huge size of these packages I'd try to link to these behemoths (well, sherpa anyway, as it has a python interface), but I'm trying to (a) use more generalized packages, (b) wean myself away from relying on the same software I've used for years and learn something new...

dd, python, google charts and amazon cloud computing

The latest RedHat magazine has a neat little article mixing dd, python, the Google Charts API and virtual machines with Amazon's version of cloud computing - the aim is a pretty pedestrian rough benchmark of disk I/O, but as a mix of computing metaphors (classic *nix command line, modern interpreted language, Web 2.0, virtualization) it is noteworthy.

Linux Journal covers number crunching with GPUs and scientific python

From a number cruncher's perspective the November 2008 issue of Linux Journal is actually quite interesting. In addition to an article on the petaflop Roadrunner supercomputer at LANL, there are two articles on number crunching with GPUs (one by Robert Farber, the other by Michael Wolfe) and one article on numpy/scipy by Joey Bernard.

Farber's article is really an advert for NVIDIA's CUDA, and unfortunately doesn't actually show any examples of doing something (this is something of a trend in LJ, as it seems to drift away from showing you how to actually code things and toward merely describing point-and-click mega-code-projects). Wolfe's article shows snippets of matrix multiplication using both CUDA and Brook, but is really more of a discussion about how to write a compiler that would automatically parallelize for GPGPU work.

Joey Bernard's article on numpy does try to scratch the surface with simple worked examples, including the use of matplotlib for plotting (or "ipython -pylab").

For a moment after reading his article I was terrified that I'd totally messed up using numpy in my python projects, as Bernard states that array multiplication in numpy (e.g. a3=a1*a2) is handled as a matrix multiplication!

However, its pretty easy to verify that unless you specifically created matrix objects (which his example did not do), then a1*a2 is an element-wise array multiplication.

To do matrix-like multiplication on numpy array objects you need to specifically do something like "a4=numpy.dot(a1,a2)" or "a5=numpy.mat(a1)*numpy.mat(a2)."

Anyway, I hope this issue is a sign that Linux Journal may get back to publishing more "hands-on" articles on general programming topics interspersed between all the Web 2.0 and database articles.

Wiki vs Google Groups/Sites

For work reasons I've been looking at wikis, as for our international collaboration will need an electronic way of outlining what we're doing, discussing it, and uploading files/images and other such info (and we're not allowed to spend money on flying furriners over here). Oh, and it has to be private and secure...

I must say I hadn't realized how horribly messy and complicated wikis are to set up, or how many non-obvious capabilities and requirements you have with the various wiki options available.

I quickly decided that I didn't want a relational database at the back end, which narrowed the field nicely, before settling on what I hoped to be a relatively simple set of wiki software: the python-powered MoinMoin. Flat file back end and powered by python - it sounded ideal to me.

Five hours later I had a full test site up, but I can't say it was as easy as I'd have liked. Nor was the security and authorization stuff at all transparent (try reading this), in fact I had no idea of whether the test site was secure or not. And wiki syntax is terrible (non intuitive!). I'm amazed wikipedia works at all.

Luckily I then discovered Google Sites ("Create websites and secure group wikis") which I must say I'm very impressed with so far. Simple Google account password protection, revision history and alerts, semi-WYSIWYG editor, comments and attachments. Maybe not industrial strength but much simpler and easier to use. Maybe we'll end up using Google Groups instead, but still much much easier than setting up a localized CMS. Good job Google!