[R-sig-Geo] Legend in density plot

[Adrian.Baddeley at csiro.au]

On 18/11/11 21:13, Albin Blaschka wrote:

> I have GPS data for three seasons of roaming livestock and to get an
> idea which areas the animals occupied and the differences between the
> seasons, I am plotting density plots for each season like this:
>
> posi <- readOGR(dsn = "PG:host=127.0.0.1 user=user dbname=mydb
> password=secret, layer = "mylayer")
>
> pposi <- ppp(posi at coords[,1],
> posi at coords[,2],xrange=c(posi at bbox[1,1],posi at bbox[1,2]),
> yrange=c(posi at bbox[2,1],posi at bbox[2,2]))
>
> densi <- density(pposi, bw = "nrd0", adjust = 0.1, kernel =
> c("epanechnikov"), edge = TRUE, diggle=TRUE)
>
> plot(densi, col = mypalette)
>
> so far, so good:
>
> BUT as there where differences between the years, the legend which is
> plotted on the left of each map shows different values/ has different
> ranges, which is bad for a comparision, so my question is:
>
> How can I change the range and the ticks in the legend of a density map?

This is a question about the 'spatstat' package.

Firstly, if you want the density maps to be strictly comparable between successive years, you should use the same smoothing bandwidth. In the code above, the bandwidth is selected automatically by the rule "nrd0". That means the bandwidth values are data-dependent and will be different in different years. You need to select a particular bandwidth value and stick to that. 

Incidentally you've used "nrd0"  which is a bandwidth selection procedure for non-spatial data. The standard methods for bandwidth selection in spatial point processes are Diggle's cross-validation method (implemented as bw.diggle) and Stoyan's rule of thumb (implemented as bw.stoyan). See the help file for density.ppp. 

Note that, because 'pposi' is an object of class "ppp", when you call density(pposi) the function that is executed is density.ppp so you should look at help(density.ppp) rather than help(density).

If X1, X2, X3, ... are the point patterns for the successive years, then I would recommend something like
        sigma <- bw.diggle(X1)
        D1 <- density(X1, sigma)
        D2 <- density(X2, sigma)
       .....

The result of density.ppp is an object of class "im". So when you plot it, you are using the plot method for "im" objects, plot.im. The help file for plot.im contains a very detailed explanation of how to control the colours and annotations. To use the same colour map for several images, you can either specify the argument 'zlim', or you can create an object of class "colourmap" specifying the colour map explicitly, and specify this as the argument "col".

Adrian Baddeley