GIS with Haskell 1

It's time to bite the bullet and do some GIS Haskelling.

My first project is to develop a simple map server in Haskell. Here are the ingredients:

• PostgresSQL + PostGIS
• Some data to put into the database. For this I sourced some Australian suburb boundaries.
• A library for manipulation GIS geometry, GEOS. In particular this provides functions to parse WKT strings from the database into geometry structures.
  
• Haskell CGI package
  
• Haskell bindings to the GEOS library.
• Extension to HaXML adding SVG combinators.

The definition of the map is specified using a set of combinators resulting in a DSL looking a lot like the MapFile format of MapServer.

The following gives a map showing the suburbs centred on the city of Melbourne.

ex0 = map (connection "host=localhost user=postgres password=postgres dbname=australia"
  `u` size 700 700
  `u` extents 144.897079467773 (-37.8575096130371) 0.16821284479996734 0.1410504416999956
  `u` layer ( table "suburbs"
          `u` geometry "the_geom"
          `u` klass ( style ( outlinecolour 255 0 0 1
                          `u` colour 100 255 100 1))
  ))

The resulting SVG file when viewed looks like:
The components that make up the map definition:

• connection - supplies the database connection parameters,
• size - is the size of the SVG output to be generated.
• extents - is the extents of the map in world coordinates.
• layer - this defines a layer. The table property defines the database table to use, the geometry property defines the column name to use. The klass parameter defines the style to use for drawing the geometry.

Layers are the key to this. Basic layers associate a geometry from the database with a style to be used when drawing the geometry.

A slightly more complex example with two layers is the following:

ex1 = map (connection "host=localhost user=postgres password=postgres dbname=australia"
      `u` size 700 700
      `u` extents 144.897079467773 (-37.8575096130371) 0.16821284479996734 0.1410504416999956
      `u` layer ( table "suburbs"
              `u` geometry "the_geom"
              `u` labelitem "name_2006"
              `u` klass ( style ( outlinecolour 255 0 0 1
                              `u` colour 100 255 100 1))
              `u` label (colour 255 255 0 1))
      `u` layer (table "suburbs"
              `u` geometry "geomunion(the_geom)"
              `u` klass ( style (  outlinecolour 0 0 0  1 `u` width 4)))
      )

This is the same as before but with a new layer that provides a thick border around the edge of all the suburb boundaries:

Next steps are to source some population data, to colour code the suburbs depending on population and to include a legend.

Financial Contracts, Haskell and Probability

This article brings together the ideas presented in the paper 'How to write a financial contract' (HWFC) and Martin Erwig's FPF module.

We are going to deal with a simple but common situation in finance - if I have a contract where I am going to receive $100 dollars in 3 years time what is that 'contract' worth to me now. How much would I pay to obtain that contract? In order to calculate the worth we need to consider what else I would do with the money and the most obvious action is to deposit it into a bank account that attracts interest.

The question is reposed then as: if I put x into a bank account then what is x if the final amount after 3 years is $100. This is easy if the interest is fixed, not so easy if it varies.

This blogpiece will provide a fragment of the implementation of HWFC that answers the above.

As this is literate Haskell some preliminaries:

> module Main where
>
> import Probability

HWFC introduces the concept of a value process which is a function from time to a random variable. We shall equate a random variable with a probability distribution and a definition of a value process is:

> type PR a = Int -> Dist a

For our interest rate model let us say that from year to the next the interest rate can either stay the same, increase by 1% or decrease by 1% all with equal likelihood. We can express this as:

> interest :: Floating a => a  ->  PR a
> interest i n  = (n *. one) i where one start = uniform [start+1/100,start,start-1/100]

The *. function allows us to repeat a random process n times. The process here is to start with an interest rate and to move to the next years rate.

If this year the rate is 10%, after a couple of years the distribution looks like:

interest 10 2
10.0    33.3%
 9.99   22.2%
10.01   22.2%
 9.98   11.1%
10.02   11.1%

Let us put that to one side and look at the contracts side of things. I will short circuit the approach in the paper and dive directly into the valuation

> data Obs a = O { evalObs :: PR a }
>
> konst k = O (\t -> certainly k)
> lift f (O pr) = O (\t -> fmap f (pr t))
> lift2 f (O pr1) (O pr2) = O (\t -> joinWith f (pr1 t) (pr2 t))
> date = O (\t -> certainly t)
>
> data Contract  = C { evalContract :: PR Float }
> cconst k = C $  \ _  -> certainly k
> when o c  = C $ disc (evalObs o) (evalContract c)
>
>
> at t = lift2 (==) date (konst t)
> zcb t x  = when (at t) x
>
> whenFirstTrue :: PR Bool -> Int
> whenFirstTrue prb = f 0 where f i = if prb i == certainly True then i else f (i+1)
>
> baseRate = 10

This is a value process such that if when the first argument is true, return the second otherwise calculated the discounted value of the first argument.

> disc :: PR Bool -> PR Float -> PR Float
> disc prb prd t = if prb t == certainly True then prd t else let s = prd t
>                                                                 t' = whenFirstTrue prb
>                                                             in discount baseRate s (t'-t)
>
> discount :: Floating a => a -> Dist a -> PR a
> discount int final time = let intspread = interest int time
>                           in joinWith (\i s -> s / (1+i/100)) intspread final
>

Lets start with a trivial example to make sure that things are working as planned

> ex1 = cconst 100

The value of this contract, as a random variable, is:

evalContract ex1 0

100.0 100.0

> ex2 = zcb 3 (cconst 100)

The value of this contact is:

evalContract ex2 0

90.90909  25.9%
90.900826  22.2%
90.91736  22.2%
90.89256  11.1%
90.92562  11.1%
90.88431   3.7%
90.93389   3.7%

The PFP library has a function to provide the expected value which can be ask of a distribution. The expected value of our contract is:

expected $ evalContract ex2 0

90.9091

Haskell

As the Haskell Wiki says

Haskell is a general purpose, purely functional programming language featuring static typing, higher order functions, polymorphism, type classes, and monadic effects. Haskell compilers are freely available for almost any computer.

Some of my contributions to the Wiki and Hackage are:

HJS - A JavaScript interpreter.
Enterprise Haskell - Requirements for the use of Haskell in the real world.
HGene - The beginnings of a geneology program in Haskell.

Lists considered harmful

A quick post inspired by the paper "Stream Fusion From Lists to Streams to Nothing at All" All programming languages include features for lists/collections. The problem with your bog standard list is that there is no tie-back to what built the list. This means that the opportunity for any optimisation that you could get by fusing the creation of the list with its use, is lost.

Conversations with a type checker

Haskell encourages a high level of thought prior to putting down characters. One feature that has been noticed is that once written a Haskell program will usually do the right thing. Haskell moves the task from punching out characters to thinking about what you are writing and, importantly, getting the types consistent across the whole program.

Haskell does not force you to specify a type for everything. This enables you to develop a function iteratively and then to ask Haskell what it infers the type of the function to be. As an example, suppose you had a higher level function that you knew the general layout of. You know that the function calls other functions but are not sure what the types of these functions are. You can get an idea of their type by writing the top level function as if the lower functions where arguments to the higher function and then asking Haskell for the type of the top level function. The type signature would include information about the lower level functions.