Spatial joins using PostGis

In a recent post, Tim talked about reproducing Spatial Joins using Hadoop. The reason why we started investigating on those topics is that we are cross querying the GBIF occurrence database with the Word Database on Protected Areas. The results of this work is already available for Spain, Madagascar and Tanzania. You can also check a particular protected area, like Sierra Norte in Spain. You can expect a following post about the widget developed in Flex and the visualization solution applied there.

In August we started the project only focus on Spain and Madagascar. After extracting the GBIF data that was less than 3 Million records and around 800 protected areas (polygons). With this "little" data I thought I could use some brute force and did not spend much time optimizing the processing of the data on PostGIS. In total it needed around 8 hours to analyze all data and place it in a DB schema that was enough for the service layer of the widget. These scripts are available here.

The project was very well received and soon we were asked to extend the analysis to all countries, all protected areas, and all GBIF data. That is around 90.000 protected areas and 130 million occurrences from GBIF. We knew from the beginning that the previous strategy would not scale to handle so much data and that the strategies for data representation would also not work with protected areas of the size of Spain and millions of occurrences inside. So this time we decided to go in two different ways at the same time: Using Hadoop to do the spatial join and adapt the PostGIS strategy to make it more scalable. This post is of course talking about the second strategy.

We wanted to work a little bit generic here and be able to support any kind of "polygon" analysis. That mean we would not be protected areas specific. So we have defined a Starting schema and a target schema databases. The starting schema is how we expect to get the data at the beginning and the target schema is the final processed schema that will be used to serve the data to the widget. GBIF gave us the primary data and WCMC the World Database on Protected Areas. The starting schema looks like:

The data from GBIF was given as a MySQL dump and imported into PostgreSQL using COPY statements. The WDPA on the other hand was given as a huge shapefile that we import using shp2pgsql. 

The main things we needed to process were stats on:

1. Find which GBIF occurrences are inside which polygons
2. Which providers have data for each polygons
3. Gridify the occurrences to be able to visualize them dynamically on Google Maps.
4. Construct a taxonomic tree for each polygon and grid.

So the target schema looks like this (click to expand):

The numbers in red are the number of records on each table after the data was processed.

This is the definition of main tables:

occurrence: A primary data record. The location of a species at a certain place in a certain moment.
site: Even if there are 130M occurrences they are only in around 3.5M distinct points. A site is a point on earth where an occurrence is. This just allows us to geoprocess only the distinct positions and save us lot of time.
grid10,grid5...: The grouping of the sites in cells at different scales, 10, 5,1,0.5 and 0.1 degree. After 0.1 we go dynamically grouping using PostGIS SnaptoGrid (another full post about it would be needed).

The whole processing is done using PostgreSQL and PostGIS with some new functions that we needed to create. The details about the whole process is described here . And the main steps are:

1) Adjust PostgreSQL: use 2GB of RAM. AutoVacuum Off. Use as many processors as possible.

2) Get the sites: No more than a grouping of the occurrences by coordinates.

3) Spatially join the sites with the polygon table: This is where it gets interesting. The fastest way we found was using the distance function and run it threaded in parallel in different parts of the site. This is a  computer intensive task so I used the 8 cores in my computer. The task takes around 4 hours. But depend s a lot on the complexity of the polygons and how good you are at distributing the load among different processors. We used a PHP script and a shell script to run several queries at the same time.

4) Delete sites that are not inside polygons and later delete the occurrences that are not in sites.

5) Create the grids: we are using our own grid system based on work previously done by Tim Robertson. The way we group the sites on the occurrences is by spatially joining the sites against tables with all the grid cells pre-created. To generate the grids Tim developed a little Java tool that create PostGIS insert statements. It is the GridBuilder class .  But if someone needs them we can provide them directly as zip SQL files.

6) Generate stats: Calculate the number of occurrences, species, specimens, observations, other basis of records, plants, animals and other kingdoms. Those stats are available per polygon, per grid, per site, per provider, per taxon, etc... Lot of queries but that run very fast.

7) Generate taxonomy for each polygon: this can be a complicate process. The idea is to generate classifications per polygon based on the occurrences that are inside these polygons. We have created a PostgreSQL function that scan the occurrence table per polygon and generate the  classification. The source code of this function is here .

And this is more or less everything. Of course lot of details are skipped here, like when to vacuum, generate indexes, etc. But you can find most of these details on the description page of this strategy .

You can find in any case all the source code for all those things on the Google Code GBIF-WDPA project. There you will find the source code of the Java SOA application, the Flex source code of the widget, the processing scripts and comments, well, everything created during this project.
we are working on a single script that will run the full process in one shot and it will be available there when ready.

Summary:
Well, after this long post my main recomendations if you need to use PostGIS to do something similar would be:
-Adjust postgreSQL settings regarding memory.
-Consider threading long operations that require a huge table scan.
-"Create table as" is much faster than updating.
-Amazon Web Services can be very slow due to I/O performance issues.
-Buy a Wii to play while you wait for results :)

Leave us a comment if you have any question or suggestions regarding this processing strategy.

And finally happy christmas everybody!

Reproducing Spatial Joins using Hadoop and EC2

Finally I got some time to continue researching the use of Hadoop for the processing of the GBIF dataset with respect to the protected areas of the world.

The problem
There is a growing number of point data, and an increasing need to perform cross referencing of the GBIF data with external sources, in this case the protected areas of the world.

The last 3 years of working on reasonably large databases (>100G) has told me that this join was not going to be an easy one, despite the obvious partitioning strategy that can be invoked here (spatial partitioning).
Therefore Javi and I set about this task by

1. Using PostGIS and the traditional approach of joining 2 tables of 150M point records with 120,000 polygon records - finding will be posted shortly
   
2. Using Hadoop and the Amazon EC2 to process the join

Using Hadoop
This is ongoing but I thought I'd post the findings so far.

Strategy 1: 

Partition by cell
Start by splitting the 150M records from one file into a file per 1 degree cell.
Loop over each polygon:

• Determine 1 degree cells covered polygon
  
• Use the input files as a merge input to hadoop and pass in the Polygon for the Map to do the contains()
  

Result  

Hadoop will not handle so many file splits no matter what configuration I try so far. Mailing list responses suggest this is not the ideal use of Hadoop and it prefers single large files, rather than many smaller ones.

Amendment

Using 10 degree cells, it does work but takes a LONG time... so long I killed it. This was expected since it is really the equivalent of only a DB join with an index on an INTEGER teg_deg column.

Strategy 2: 

1-Degree in-memory index of polygons
So this approach goes through the join the opposite way. The polygons are loaded into memory and a

HashMap<Integer, Collection<Polygon>>

for each Map operation to use.  

The logic is effectively give me polygons intersecting the 1 degree square for this point, and then test them properly. This worked ok, but requires a LOT of memory for the polygons to be held.

Result

About 6 hours to process 60,000 polygons.

Strategy 3: 

RTree index of polygons
Very similar to strategy 2, this approach indexes all the polygons by RTree. Thus for each point, the index provides only candidate polygons whose bounding box contains the point, and the number of candidate polygons to test is far smaller than that of strategy 2.

Result

45 minutes for a 130M x 30,000 cross reference but requires the larger EC2 instance sizes. The RTree index means it will need processed in batches.
Hadoop requires tuning, just like a database does for large joins. Here is the final tuning for strategy 3:

• mapred.child.java.opts=-Xmx1G // for the R-Tree index
• mapred.output.compress=false // for easy reading of output
• mapred.tasktracker.map.tasks.maximum=4 // we have 7G per node, so 4G for 4 Map tasks is ok
• mapred.job.reuse.jvm.num.tasks=-1 // Reuse the Map as it has the index built (takes 90 seconds for startup)
• dfs.block.size=134217728 // use larger blocks - in retropsect this was a probable mistake and will be removed for next test

More to follow...

Thanks to the Hadoop-core Users mailing list for their insights, and Andrea of OpenGEO for pointing me at the JTS RTree implementation

Do you know any species distribution project?

I am trying to create a general view on how species distributions are handle in different projects world wide. I am specially interested in project where they have maps and are published on Internet. So I would like to ask you for links to projects related to Species Distributions.

I am also interested in different ways people handle occurrence status or absence. 

And finally I also would like links to any mapping application that you have seen related to Species distributions. And if possible comment what you liked.

Please post your links as comments to this post or send them to me (jatorre@gmail.com). I will post in the future the results from this little research.

Thanks!

Flickr creating polygons out of points, what about species distributions?

This weekend I have found a very interesting post from the Flickr programmers called The shape of Alpha. What this people are basically doing is generating shapefiles from the geotagged pictures in Flickr. When a user geotag a picture, provide lat/long, they reverse-geocode it to find the neighborhood, the town, the country etc where this coordinates are.  They store this information in the database (as  Where On Earth(WOE) IDs) and this allow for much better performance when later people search for pictures in "Madrid".

But this post talks about a different idea, and it is awesome. Now that they have all this geotagged pictures and they know they are in a particular neighborhood or country (WOEID) they are creating shapefiles from these areas based on the coordinates of the different pictures! Maybe is easier to explain it with a picture:

The one in the left is a London polygon and the right for United States. They were created by aggregating all the coordinates of the different pictures found in London and in the US.

Again. They take all the coordinates, POINTs, in the DB for a certain WOEID and with them they generate a polygon. If they have enough data it looks more or less ok. They have a good discussion about the quality of the polygons and the threshold set to generate them. I like very much the introduction they make to  Alpha shapes and the links they provide.

Additionally they also provide the source code of the tool they use to process the data. It is called Clustr.

Now, lets apply this concept to Biodiversity, if you haven't already figure it out. Think of their geotagged pictures as PRIMARY DATA, the WOE as SCIENTIFIC NAMES and the polygons you get out of them as DERIVED SPECIES DISTRIBUTION polygons.

I have been investigating this for a long time already, specially trough the Biodiversity Atlas project that is now stopped, but hopefully will start soon. We can take their source code and apply it to GBIF data and generate "derived, unchecked and uncompleted" species distributions based on GBIF data on a massive way! And the same they try in Flickr, the more primary data we get into the system the better the distributions will start looking like.

But there is another idea... why dont ask Flickr to not only process their polygons based on WOE Ids but also on tags? So if we tag pictures in Flickr with scientific names, or better GUIDs, they can then try to generate by themselves the distributions.

I particularly dont see Flickr as the best place to handle species distributions discussions, and will, hopefully, try to convince at least one big biodiversity project to let me try this way. Most of you can probably imagine the incredible API we can create once we have a lot of species distributions accessible in such a system. I will write another post about it soon, but think of:

• Which species could live in my garden?
• Which species habitats this new road will cross? 
• Where should I create a new Protected Area to preserve as much biodiversity as possible?
• What species could I find in the track I will hike this weekend?

So maybe instead of so much niche modeling projects we should start thinking on how to manage the vast amount of primary data (and other sources) we have already  and how to curate it and complete it. I dream of a scientific community joined together to create a complete information system to know where species are. Imagine like a Wikipedia but for Species Distributions.

Word clouds for TDWG 2008

Some days ago I found a great application called Wordle. The application let you generate "word clouds" from text that you provide. It is very similar to what a tag cloud is but with the difference that it analyze the words in a document and not just the tags. The application analyze the frequency of word occurrences.

I gave it a try and I thought it could be interesting for others. For example here is the Word Cloud for the Proceedings of TDWG 2008. Basically I just copied all the text from the PDF and gave it to Wordle for analyzing.

Click on the image to view it in their Java applet.

Then I did the same but with Markus post about TDWG 2008. And here is the result:

They look similar but Markus I think has a higher geek level :D

Finally I could not resist and checked how it would look this blog and here it is:

Seems that the words that use most are Google and Data, even more than biodiversity itself (10 counts for Google and 6 for biodiversity). 

Finally, if you are interested in displaying word clouds in maps there is a recent post for one of my favorites blogs, cartogrammar. I dont see much use on it apart of that they look nice, and they indeed look nice!.

TDWG 2008 in Fremantle

The great TDWG conference has just finished this year in Fremantle, Western Australia. Me and 180 other biodiversity informatics people have attended and I thought I would point out the most interesting things from my side.

Darwin Core

The most exciting thing to me is the chance to come up with a new Darwin Core soon for ratification that resembles more Dublin Core and consists of possibly 3 namespaces.

While working on the IPT it seemed that most information we are dealing with in biodiversity informatics is centered around 3 entities only: Taxon, Occurrence and SamplingSite. It would make sense in my eyes to separate elements of a new DarwinCore standard according to those core entities. The normative standard will be decoupled from the implementation technology and only consist of natural language definitions with URIs for the elements of the standard. Specific application schemas (called profiles in DC) will then define/recommend how to use DarwinCore within an XML, RDF, XHTML, OGC application schema or tab file environment. Even datatyping can be left to the respective schemas and e.g. dates can be expressed in their native formats. I will give examples of the different representations soon in another blog entry.

Combining these 3 core entities with the notion of one-to-many "extensions", the star schema of the IPT, one can handle quite a rich definition of data. Extensions for multiple identifications can be added to the occurrences, SPM like species descriptions, geographic distribution or species status for invasiveness to the taxon. And still the whole standard and exchanged data can be extremely simple! Btw, simplicity was probably the most mentioned idea in TDWG this year (a really nice talk was about creating a Species Index with SiteMaps by Roger Hyam in our Wild Ideas! session). Controlled vocabularies like BasisOfRecord or Ranks should be expressed in simple ASCII files like the ISO country code one, with a code, label, definition and examples. This also allows for easy translation in different languages because the codes stay the same.

Integrated Publishing Toolkit

Surprisingly many people asked for IPT demos after screenshots were shown in some GBIF talks. Me and Tim did demos nearly every day and the publishing tool was well received in general. Alpha testers (mainly usability) for the public instance at GBIF were gathered and new ideas arose or vague ideas materialised, e.g. validation/annotations:

 * validation can be done through external services that adhere to a simple API. The validaion would be asynchroneous, sending a token, link to the full dataset (either dwc xml or tab file dumps) and a callback handler. Once the validation is done the handler would receive the token and a link to the validation report, an xml file that contains annotations (unstructured text) about records together with some probability and potential suggestions (list) of property changes.

 * provide API to push datasets into the IPT via REST service

 * BCI collection and institution code validation and lookup of GUIDs during mapping.

TAG & GUIDs

Greg Whitbread will likely be leading the Technical Architecture Group. Refining and shrinking the core ontology (owl) was seen as an important issue over the next years.

Originally I had planned to question the uptake of LSIDs again, being in favour of PURLs since the beginning. After long debates that some of us seemed to have experienced before, we came to the following conclusions though keeping up the LSID recommendation:

 * LSID look much more stable in printed publications. That means they should really be resolvable through proper LSID resolution via DSN SRV records

 * LSID are on the agenda of many projects already

 * Proxying LSID with http removes many troubles, especially when used with RDF. The strong recommendation is to always use the proxied version in RDF abouts.

 * Changing the domain used within the LSID might cause problems (e.g. name change of an institutions). It might be better to have a single central LSID authority

 * pure UUIDs with a central resolver would be great. The resolver would have to know about all existing UUIDs in our domain, but it could be mirrored and sharded easily. Maybe something to think about

 * central PURLs and LSID (only) require a registration of services/authorities

Names

Great introductory talks were given by Rich Pyle about zoological and botanical nomenclatural differences and by Nico Cellinese about the phylocode. If you always felt to be on shaky grounds with nomenclature or taxon concepts you should really watch these talks!

Identifying good images on Google cache for scientific names

I have been working on ways to represent taxonomic trees in a more intuitive way for non-biologist. The best I have found until now is to provide a contextual image together with the name. This helps a lot specially on higher ranks where most of the names are very unfamiliar for most of people, like me for example.
The problem is from where to get images about all scientific names. Well the best I managed to find is Google AJAX API. This is the content you get when you search for images on Google. The service is fast and reliable. But, it has one problem, it is not content aware. Well, they are trying, check this and much better this, but still it is not. So sometimes the results you get back just dont make any sense. There are some very pornographic examples of this, but to keep it kids friendly check the Phylum Labyrinthulomycota. The first image you will get is from this nice researcher:

Well, that is not of great help to get an idea of what this phylum is about. Therefore I had been for a while thinking in doing a little application where people can help me to select pictures that actually help to get an idea of what is behind a cold scientific name. And today I had a little bit of time and wanted to deploy something on AppEngine.

So here I am presenting a very simple application to ask for collaboration on this task. There are only 13 million names that I need to find a picture for, but I think I have enough friends :D

http://snimagefinder.appspot.com/

The application is very simple and I have not added any visual effects, but I just wanted to give it a try.

The rankings I get from this will be released soon in an API for anybody else.

Things I would like to incorporate:

1) Make it a game, more precisely a GWAP "game with a purpose". Like Google Image Labeler or any other toy from Luis von Ahn. (I wanted to link to his site for sooo much time).

2) Allow Multiple rankings. Now names only get evaluated once so a malicious contributor can ruin all this.

3) More sources, specially Flickr.

4) Nice UI so that Sergio, our new blogger here, is happy :)

When I get a decent amount of reviews from people I will post some stats and an example application.

Come on send the link to everybody and help everybody to understand better scientific names!

http://snimagefinder.appspot.com/