[R-sig-Geo] Plotting Kernel results with ggplot2

Diego Pavon diego.pavonjordan at gmail.com
Mon Feb 29 16:01:05 CET 2016


Hello Barry,

Yes. What I meant was to have the contours that have 50% and 75% of the
points. But as far as I understood, this is what a kernel analysis does,
right? Gives you the minimum area (smallest polygon) containing the level
specified (50% of the data points, 75% of data points...)?

If I do the kernel analysis with the kernelUD() function in the adhabitatHR
package, I can easily get the contour and the polygon to be plotted even on
a map from 'library(rworldmap)':

xy <- data.frame(x = runif(30), y = runif(30))
coordinates(xy) <- ~ x + y
plot(xy)
class(xy)

library(adehabitatHR)
# calculates the Kernel
kernel.sim<-kernelUD(xy, h = "href", grid = 100,
                          same4all = FALSE,
                          kern = c("bivnorm", "epa"), extent = 1,
                          boundary = NULL)

image(kernel.sim)
str(kernel.sim)
summary(kernel.sim)

# gets the polygon of the kernel to be drawn
ver.sim50 <- getverticeshr(kernel.sim, 50) ## 50 %
ver.sim50
ver.sim75 <- getverticeshr(kernel.sim, 75) ## 75 %
ver.sim75

plot(ver.sim50, add = TRUE, lty = 2, lwd = 3)
plot(ver.sim75, add = TRUE, lty = 3, lwd = 3)
points(xy, pch=16, cex=0.5, col='royalblue4')

# GET THE COORDIMATES OF THE CENTROIDS (center of the polygon)
coordinates(ver.sim)

#Draw the centroids
points(coordinates(ver.sim50), pch=16, cex=1.5, col='olivedrab4')


However, I can't plot these contours, kernels, polygons, etc in ggplot2. Of
course ot would be handier if there would be a way to just plot in ggplot2
the results from the kernel analysis shown above.

Do you have any idea about this?

Thank you very much

Diego







2016-02-29 12:21 GMT+02:00 Barry Rowlingson <b.rowlingson at lancaster.ac.uk>:

> I think you will have to take off the ggplot2 training wheels and do
> it another way. Here's the outline:
>
> 1. Use kde2d to compute your kernel on a grid - you'll have to choose
> the bandwidth and grid size, ggplot2 usually makes those decisions for
> you. Something like:
>
> k = kde2d(xy.df$x, xy.df$y,h=0.3,n=200)
>
>
>  2. Convert the gridded density to contour lines with
>
> cc = contourLines(k,levels=c(.25,.75))
>
>  with whatever levels you want. Now I'm not sure what you mean by
> "50%" and "75%" - do you mean contours at that fraction of the peak
> value? Find the max of k$z if you do and work with that. The absolute
> values of k$z are dependent on the cell size since it should integrate
> to 1 over the area. Or do you want contours that contain 50% and 75%
> of the points? That's not quite so simple.
>
>  3. Convert the `cc` output to something ggplot can handle.
> `contourLines` returns a list with one element per contour line (loop
> or segment) with the level, x, and y coordinates. Suspect you need to
> `cbind` all the coordinates with the level and the index to get
> something like:
>
>  index level x y
>   1      0.75 0.123 0.545
>  [etc for contour 1]
>   2     0.5 .5345 .9123
>   [etc for contour 2]
>
> and then you can feed this to geom_path, grouped by index and coloured
> by level.... There's probably a solution for this already somewhere...
>
> Barry
>
>
>
> On Mon, Feb 29, 2016 at 7:34 AM, Diego Pavon
> <diego.pavonjordan at gmail.com> wrote:
> > Hello list
> >
> > I am trying to do a Kernel analysis of bird ring recoveries and I would
> > like to plot the kernel contour in ggplot. I have found that the command
> > 'stat_density_2d()' can do this kernels. However, when I plot it, it
> draws
> > all the kernels from 10%... 90%.  See example below:
> >
> > # Simulate data
> > xy <- data.frame(x = runif(30), y = runif(30))
> > coordinates(xy) <- ~ x + y
> > plot(xy)
> > class(xy)
> >
> > # Transform SpatialPoints into a data.frame to be used in ggplot2
> > xy.df <- as.data.frame(xy)
> > head(xy.df)
> > str(xy.df)
> >
> > # Make the plot
> > p <- ggplot(data=xy.df, aes(background="white"))
> > p <- p + geom_point(aes(x=x, y=y),
> >                     data = xy.df,
> >                     size = 4)
> > p <- p + xlab("X") + ylab("Y")
> > p <- p + theme(text = element_text(size = 15))
> > p
> >
> > # This makes the kernels
> > p2 <- p + stat_density_2d(data = xy.df, aes(x=x, y=y), col = 'red',
> alpha =
> > 0.6)
> > p2
> >
> >
> > However, I am only interested in drawing the lines for the kernel 50% and
> > 75% and not all of them.
> >
> > Does anyone has an idea how to select only those lines (in my case 50 and
> > 75%) to be plotted and omit the rest?
> >
> > Thank you very much in advance for you time and help.
> >
> > Best
> >
> > Diego
> >
> >
> > --
> > *Diego Pavón Jordán*
> >
> > *Finnish Museum of Natural History*
> > *PO BOX 17 *
> >
> > *Helsinki. Finland*
> >
> >
> >
> > *0445061210
> https://tuhat.halvi.helsinki.fi/portal/en/persons/diego-pavon-jordan%288d5db37c-eddd-4fca-92cd-9c9956a42b4a%29.html
> > <
> https://tuhat.halvi.helsinki.fi/portal/en/persons/diego-pavon-jordan%288d5db37c-eddd-4fca-92cd-9c9956a42b4a%29.html
> >
> http://www.linkedin.com/profile/view?id=170617924&trk=nav_responsive_tab_profile
> > <
> http://www.linkedin.com/profile/view?id=170617924&trk=nav_responsive_tab_profile
> >https://helsinki.academia.edu/DiegoPavon
> > <https://helsinki.academia.edu/DiegoPavon>*
> >
> >         [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > R-sig-Geo at r-project.org
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>



-- 
*Diego Pavón Jordán*

*Finnish Museum of Natural History*
*PO BOX 17 *

*Helsinki. Finland*



*0445061210https://tuhat.halvi.helsinki.fi/portal/en/persons/diego-pavon-jordan%288d5db37c-eddd-4fca-92cd-9c9956a42b4a%29.html
<https://tuhat.halvi.helsinki.fi/portal/en/persons/diego-pavon-jordan%288d5db37c-eddd-4fca-92cd-9c9956a42b4a%29.html>http://www.linkedin.com/profile/view?id=170617924&trk=nav_responsive_tab_profile
<http://www.linkedin.com/profile/view?id=170617924&trk=nav_responsive_tab_profile>https://helsinki.academia.edu/DiegoPavon
<https://helsinki.academia.edu/DiegoPavon>*

	[[alternative HTML version deleted]]



More information about the R-sig-Geo mailing list