A Middle-aged Geek's Guide To The Arkestra : First Movement

Intro

"Arkestra is an intelligent semantic Django-based CMS for organisations and institutions". Which is seemingly exactly what we want for our department web pages.

The showcase site, Cardiff University's School Of Medicine is rather impressive.

Arkestra is built on Django-CMS, and is open source. So I decided to play around with it. This blog entry will relate my progress. All this work is being done on an Ubuntu 10.04 box.

Installation Packages

First, create a new virtualenv and get into it. Install the required packages using pip, and the cutting-edge stuff from source. Some of these packages are necessary for django-cms and Arkestra, and some are just necessary for the Arkestra example. I'm not sure which is which, so I grab everything.

virtualenv --no-site-packages Test
. Test/bin/activate
cd Test
pip install PIL
pip install django
pip install BeautifulSoup
pip install django-typogrify
pip install pyquery
pip install -e git+http://github.com/stefanfoulis/django-widgetry#egg=django-widgetry
pip install django-filer # gets easy-thumbnails, django-mptt

pip uninstall -y django-mptt

pip install django-cms
pip install django-polymorphic
pip install django-appmedia
pip install -e hg+https://bitbucket.org/spookylukey/semanticeditor#egg=semanticeditor
pip install -e git+git://github.com/ojii/django-classy-tags.git#egg=classytags
git clone git://github.com/evildmp/Arkestra.git src/arkestra

I tried to do all this with a pip requirements file but failed - the problem seems to be that django-filer will pull in django-mptt, but django-cms can't work with it.

For reference, here's the output of pip freeze telling you what versions of what I have this working with:

BeautifulSoup==3.2.0
Django==1.3
PIL==1.1.7
South==0.7.3
distribute==0.6.10
django-appmedia==1.0.1
-e git://github.com/ojii/django-classy-tags.git@8ad613b3bc0310ba9afb206136cfadbdfa8e6b01#egg=django_classy_tags-dev
django-cms==2.1.3
django-filer==0.8.2
django-polymorphic==0.2
django-typogrify==1.2.2
-e git+http://github.com/stefanfoulis/django-widgetry@57bbf529509e805a74b3e8f36e6938bef591e505#egg=django_widgetry-dev
easy-thumbnails==1.0-alpha-17
lxml==2.3
parti-all==0.0.7.22
pyquery==1.0
-e hg+https://bitbucket.org/spookylukey/semanticeditor@3a40b4766f32262be2cfe131cf20a19703bfc339#egg=semanticeditor-dev
smartypants==1.6.0.3
textile==2.1.5
wsgiref==0.1.2

Config and Run

Head into the arkestra directory that git got for you, and try running the server:

cd src/arkestra/example
export PYTHONPATH=.. # to find the arkestra modules
python manage.py runserver

It should run, but you will fail to log in because there's no user accounts. Also, the example.db file isn't compatible with the version of django-cms and possibly other packages. So you can move it and run syncdb again:

mv example.db example.db.orig
python manage.py syncdb   # create db tables, create admin user
python manage.py runserver

Now you can log in and create Page objects, and view them on the site. Don't forget to tick the published and navigation boxes if needed. To add content, you add a Text/Layout plugin and put the content in there.

Some 404s

Watch your log output. At some point you may see some 404s caused by some magic autocomplete jQuery javascripts not being found. Try this from the example/ folder:

cd media
ln -s ../../../arkestra/arkestra_utilities/static/ javascript

That should stop the 404s.

Loading the Demo Site

Note that there's a JSON-encoded version of the database supplied in example.json. In theory all you need to do is:

python manage.py loaddata example.json

over an empty database and get it all going. But there's those pesky database differences to deal with. I did go through it editing the changes and got it mostly working, but you'll probably want to build your Arkestra site from scratch. Which is what the second movement in this piece will probably be about.

Thoughts on OSGIS 2011

OSGIS is the UK Open-Source GIS conference. This is the third year it has run, and the first time I've made it there. Nottingham University has been its home every time, mostly due to the presence of the Centre for Geospatial Science on campus, and the enthusiasm of Suchith Anand and Jeremy Morley, who seem to be running on free energy all the time.

The meeting follows the standard template for scientific and technical conferences - there's some workshops, some invited keynote papers, and some general presentations. There's no poster session, and no novel ideas like lightning talks or group discussions - except perhaps in the Open Gov Workshop, but I wasn't there. Not that group discussions don't happen elsewhere, but only over coffee, great lunches, and beers.

Links to the agenda and presentations are all on the Web Site so I won't go into details here. The talks were also webcasted live and recorded, and should be available soon.

The workshop sessions on QGIS were run by Faunalia who let us take away their 8Gb memory sticks packed with a bootable Linux. Very handy. There were also lots of OSGEO Live DVD discs kicking around too, used for the workshop on Ordnance Survey data. Clearly open-sourcers are finding more ways of working in labs full of Windows boxes. I once ran a workshop where we used Virtual Box and an OSGEO Live DVD image to hide the horrors of Windows XP.

So what was hot or not at the conference? Here's my impressions:

On the desktop:

• Hot hot hot: Lots of mentions of QGIS.
• Also Hot: gvSIG (a day of workshops), GRASS here and there, OpenJump poked its nose in.
• Not Hot: No-show GIS packages: uDig, SAGA.
• Hot: OpenLayers and GeoExt.
• Not hot: The Cloud. Its not the answer to everything.

On the server:

• Hot: PostGIS continues its dominance.
• Not Hot: Mentions of Oracle Spatial were usually accompanied by apologies. Usual excuse was 'we had the license anyway'.
• Not Hot: No sign of NoSQL DB technology - where's the MongoDB Spatial people?
• Hot: WFS, WMS, WCS and even WPS becoming almost commonplace.

In the hands of the programmer:

• Warm: Only one statistics package mentioned - R of course - although the Orange Toolbox gets a look-in.
• Hot Programming languages: Python, Java, JavaScript, C, R, SQL. I think I saw some Perl.
• Cold Programming languages: C#, Ruby, Visual Basic/C.
• Warm and fuzzy feelings from: Drupal, Django, ExtJS.
• Cold and dead to us: Sharepoint. IIS. .net. Silverlight. Flash.

Hardware:

• Hot: Android phones - nearly everywhere!
• Cold: iPhones - hard to spot. Either nobody had them or were ashamed to show them.
• Cold: Tablets. Think I saw one Samsung Galaxy Tab but nobody was showing off iPads. Or iPhone Maxis as they are sometimes known.

No great surprises there.

A few people were tweeting useful things during the conference, and so I should apologise to Jeremy Morley for not saying hello IRL and only communicating via twitter. I will try and submit a paper on doing GIS in R for next year. It will be more than 140 characters.

Also great to finally meet Suchith - the force behind lots of interesting things going on at Nottingham.

Tyler Mitchell rolled in from across the pond and did well against the jet lag to present a keynote speech and network like crazy for two days. OSGEO are very busy both globally and locally.

In all the conference did feel like a mini FOSS4G, but with a local feel - but not a restriction that cut down the quality of the talks or the participants or the coffee chat or the lunch discussions or the beer conversations. See you next year!

Reading Layers from a Qgis Project File.

Some notes on reading individual layers out of a Qgis Project file. This should end up in a 'import layers from a project file' plugin.

1. get the xml into some text:

xml = file("project.qgs").read()


2. make a blank QDomDocument and set the content:

d = PyQt4.QtXml.QDomDocument()
d.setContent(xml)


3. Now find the maplayers:

maps = d.elementsByTagName("maplayer")


4. How many?

maps.length()


5. Now read the fourth layer (0-index) into running Qgis:

QgsProject.instance().read(maps.item(3))

Wrap this in a nice gui that tells the user what the map layers are, lets them check boxes, import the layers, marks the current project as dirty. Jobzagoodun.

OpenLayers 2.10 Beginners Guide Book Review

A man is out driving in unfamiliar countryside. Inevitably he gets lost, so he pulls over to ask a local farmer for directions. "Oh no", says the farmer, "I wouldn't start from here".

It's an old joke (and not a very funny one), but one that springs to mind when people ask certain questions on mailing lists. Someone with very little computing experience will appear on the OpenLayers mailing list (or IRC) and say "How do I put a map on my web page?". My answer would be "Well, I wouldn't start from here".

First, I'd say, learn programming using a nice friendly language - python perhaps. Then get a basic understanding of how client-server systems work, and how the web is implemented via the HTTP protocol. Then follow the 20-something year evolution of HTML - from its humble beginnings through to its happy marriage with CSS and their successful menage-a-trois together with Javascript. Ignore its affair with Java - that was a moment of crazed passion, and Java is now out of the HTML/browser party and busily running the "enterprisey" affairs of corporations.

With a knowledge of HTTP, HTML, CSS and Javascript, you are ready to put maps on your web pages. And I'd choose OpenLayers to do it.

Erik Hazzard's book, "OpenLayers 2.10 Beginner's Guide", published by Packt Publishing, is designed to get people with very little computing experience all the way through that process to the point of quite complex map production. He has to cover everything, and does so pretty well (although I'm not sure about his exposition on Object-oriented programming). Although he doesn't teach programming, he encourages users to cut-and-paste code and then experiment.

The sample chapter will show you how the book works - simple explanations, example code, and open-ended questions and mini-assignments. The Firebug debugger is introduced early on in the book, and I think this is a great idea. Firebug didn't take long to repay me the time I took to grok what it could do. Erik sprinkles it liberally on the examples, since interactive inspection of a program is a great way to really see what it's doing.

Inevitably with a book like this, large chunks of its 372 pages will be similar to dumps of the API documentation from the OpenLayers web site. Erik has fleshed out a lot of the terse API information to include fuller examples to make the book more than just that - I can imagine hardcopies of this book with a multitude of flapping yellow sticky notes for their owners to get rapid access to useful pages. Packt Publishing will sell you the PDF version for those who like to search, and the paper version for those who like to stick.

By the end of the book you should have the knowledge to understand how the sample Flickr map application developed in the last chapter works, and you'll also have a handy reference to a large percentage of the features of OpenLayers. The one glaring omission is that there is no mention of the OpenLayers "Popup" class, typically used for "info-balloons" when clicking on map features. But if you can work through the book, then it should be no problem to work out how to use them from the API documentation.

There seem to be very few books on OpenLayers on the market at the moment, so if you want to make maps and would like a high-quality step-by-step guide to it, then this is the book for you.

Note that OpenLayers, like all web technologies is evolving fast, so watch out for updates and get the latest version of the book. A recent OpenLayers hackfest concentrated on developing functionality for mobile devices, so expect to see more OpenLayers applications for your smartphone soon. With Google insisting on presenting advertising along with their map images, it might be a good time to switch to OpenStreetMap for your map background and OpenLayers for your mapping client.

Barry


Note: I was provided with a PDF copy of the book for review, and I have no connection with the author or publisher.

Tolerance and bridges

So I'm modelling the H1N1 Swine Flu outbreak in Glasgow. One of the questions is how different regions (specifically 'datazones') within the study area interact with each other. You'd expect zones next to each other to interact with each other more than those on opposite sides of the city, for example. So our first idea is to use the distance between centroids as an explanatory variable and a distance-decay parameterisation for the effect. That works quite nicely.

But the real world isn't quite like that. Glasgow has a river running through it, splitting our study region into north and south parts. I wanted a new distance metric to reflect this partition.

So I came up with this idea - that the distance between two datazones is the minimum distance between the centroids of the graph formed by adjacent datazones.

First problem of that is computing the adjacencies between datazones. In R, I can just use poly2nb from the spdep package and get that. But I found that the graph created wasn't connected over the study region. Further investigation showed that polygons that clearly looked adjacent weren't being considered as adjacent. Here's a screenshot that shows three regions that look adjacent:

Now why wasn't poly2nb saying they were connected? I zoomed in for the answer:

Check the scale bar. These polygons are about a hundredth of a millimetre apart. That's a very tiny no-man's land. More like a "no-fly zone".

I figured I needed some kind of buffer on the polygons, and started looking into the rgeos package, recently released onto CRAN. It'll do that, and all sorts of things, but then when I tried working out adjacencies using its gIntersects function I discovered it worked to those tolerances already. By some magic the apparently adjacent regions were adjacent again.

However the graph was still not fully connected - it was composed of two subgraphs - the river is wide enough, and runs fully from east to west, that no fuzzy tolerances could really overlap the polygons. And that would probably be the wrong thing to do anyway.

So I loaded the regions into Quantum GIS and then loaded some Google Maps as underlays using the marvellous Google Maps plugin. I created a new vector layer of lines and started adding extra connections where bridges cross the Clyde:

These lines, shown in red above, were saved as a shapefile, then loaded into R and their start and endpoints overlaid on the datazones. That gave me some extra edges to add to the adjacency graph, all the graph stuff being handled by R's igraph package (including the distance-weighted minimum distance algorithm along adjacent datazones).

So now I have a completely-connected graph of datazones, which I can convert into a full distance matrix that reflects the geometry of the study region. The university's High-End Cluster is currently churning out the number for the model fit, and the results will hopefully be published in a journal sometime...

Thanks to: R, Quantum GIS, Google Maps plugin authors Ivan Mincik and Nungyao Lin, R package developers for sp, spdep, rgeos, igraph.

Mapping The Kendal Mint Quake

About 11pm last night my wall of photographs rattled. It's not unusual, it seems to happen whenever anyone is having a washing machine delivered and a truck crawls past my house. But a minute after that I was on the computer and someone local was tweeting about her house shaking. Earthquake? A few twitter searches revealed twitterers in Cumbria and North Lancashire had indeed experienced a small earthquake.

So I did what all statisticians would do - go looking for data. Tweets can be geo-coded, but getting reliable info might be hard. I checked the USGS earthquake site, but it's not real-time. But then I found the British Geological Society earthquakes page, with true real-time reports from a bunch of sensors around the country. Under the Monitoring menu there's the Real-Time Data option.

Here you can choose the station and a day and get a seismogram. Here's what the one in Keswick (KESW) was showing:

BIG spike. It went off the scale of the plot. I looked at a few other plots and discovered most of them had recorded it too. And I noticed the quake had hit at slightly different times. Maybe I could use that information to work out the location of the quake.

At this point the BGS site was getting slow as worried Cumbrians checked to see if it was an earthquake or if the local nuclear power station had exploded. Anyway, it was bedtime, so I left it to the morning, as something to do in my lunch hour.

I couldn't find a way to get the raw data, so I decided I'd get what I could from the images. A quick script pulled all the GIF images from the sites, except for a couple that seemed to be off-line. Another script pulled all the station information book HTML files and from there I got the locations.

I then stepped through each of the 24 images, identifying the point near 11pm where the tremors started, and entered the data.

Now I had the location and time of the points where the quake had been sensed. This was all converted to OSGB coordinates using rgdal in R.

I didn't have a zero-time for the quake or a location. A quick plot did show the quake had to be somewhere in Cumbria.

I don't know how professional seismologists work out earthquake epicentres, so I came up with something quick and dirty.

Assume the velocity of the quake waves are constant across the UK. Then if the quake was at location (x,y) then a plot of distance against time should be a straight line, with all the points on it. Some quick experiments showed I could fit a straight line then use the sum of the residuals squared as a measure to be minimised over (x,y). A few lines of code later and I had a coordinate!

> fit$par
[1] 307499.9 506035.7

By now the quake had hit the news and the USGS web site was showing it:

How did that compare to my point? I converted it to OSGB coordinates:

> coordinates(usgsos)
     coords.x1 coords.x2
[1,]  328994.6  500054.7

20km further east, and 6km south. Not bad considering I had to eyeball the time from the seismograms and various other assumptions. I don't know what kind of uncertainty the USGS put on their predictions, and I can't put an uncertainty on my location either since there's not a proper statistical model behind this. I can give you the data if you have a better idea - or you can scrape it again from the web site! Here's a comparison of the epicentre predictions on a map:

Here's what the distance/time plot looks like for my predicted epicentre:

And from the slope of this we can get the speed of travel of the quake waves - about 7,500 km/s - which is about what a quick google search on earthquakes tells me.

Footnote: I think on the bike home tonight I worked out a way to do this statistically. Suppose the earthquake is at (x,y), at time t0, and the waves move with speed v. Then suppose the time of arrival of the quake at a point at distance d is Normally distributed with mean time t0+d/v and variance k(d;theta), where k() is some increasing function of d. This allows us to consider the increasing noise in the system for points further away. Fit the whole thing by maximum likelihood and we get estimates (with errors) for location, time of origin, and wave speed. Something for next lunchtime...

Speeding up some R Code

Part of my job involves looking at code written by research students and trying to make it do other things. Often a student is focused on their application area or particular dataset, and so the code and data are inextricably linked. Add another observation to the dataset and you'll find yourself having to change 3467 to 3468 in a dozen places. And maybe find some 6934s and change them to 6936s. The other problem is that research students are mostly new to serious programming, and so still need to develop both intuition and formal skills.

Today I've been working on such code for predicting epidemic outcomes. It was taking 14 seconds to run, which isn't too bad (especially considering some other work I'm doing takes 6 hours using OpenMP on an 8-core machine), but in my inspection I'd spotted some odd R code.

exp(var%*%t(coef[1,]))

The oddity here is that var is one-dimensional and so is coef[1,]. This matrix multiply is really just something like exp(sum(var*coef[,1])). So I tried that, and the code took 6 times longer to run. Matrix multiplications are efficient in R, even in one dimension.

So I switched back to the original code. Time to get profiling. This is easy in R.


Rprof("run.out")
dostuff()
Rprof("")
summaryProfile("run.out")


This told me most of the time was spent in the t function, within the function I was playing with earlier. Then I noticed that pretty much wherever coef was used, it was transposed with the t function.

I refactored the code to do t(coef) once, and then pass the transposed matrix around. That sped the code up from the original 14 seconds to 3 seconds. Sweet.

Note that at each stage I check that the new output is the same as the old output - it can be useful to write some little test functions at this point, but doing formal test-driven development can be a bit tricky for statistical applications in R.

There's still a bit of work to do to make this code a bit more general purpose. Currently it works on predictions for one particular country - we want to do predictions over half of Africa!