Seite 1 von 1

Save/Export specific contour levels for a KDE (package 'ks' in R) as ESRI shapefile

Verfasst: Di Sep 03, 2019 2:30 pm
von raoulreding
Hallo, da ich das Ganze bereits umständlich in englisch verfasst habe, hoffe ich es ist OK, den englischen Test einfach einzufügen. Die Antworten können gerne auf deutsch oder englisch sein. Danke!

Hello! I'm performing a kernel (KDE) home range analysis of an individual with the package 'ks'. I'd like to avoid using Geospatial Modelling Environment (GME) for implementing the output in ArcGIS, because I'm using QGIS which is not compatible with GME. I know that 'adehabitat' is much more easier to handle, but it doesn't supply the plug-in bandwith selection, which is essential for the species I'm researching. GME only works with ArcGIS and R and uses the 'ks'-package. It calculates exactly what I want, so there has to be a code in the background, which makes it possible to do it manually in R only. Here's what I want to calculate:

Here's the data:
gps_data.png
gps_data.png (9.37 KiB) 2613 mal betrachtet
And here's the kde-code:

> gps_hpi<-Hpi(x=gps,pilot="samse",pre="scale")
> gps_kde<-kde(x=gps,h=gps_hpi,compute.cont=TRUE)
> contourLevels(gps_kde,cont=c(50,75,95))
> contourSizes(gps_kde,cont=c(50,75,95),approx=TRUE)
> plot(gps_kde, display="filled.contour", cont=c(50,75,95))
> legend(460000,5660000,c('50 %','75 %','95 %'),fill=c("red","orange","yellow"))
gps_plot.kde.png
gps_plot.kde.png (5.99 KiB) 2613 mal betrachtet
These contours are exactly what I want to export as a shapefile, but there is no reasonable possibility described in the internet to do this. I've already tried many different ways to do it, I also reached to export some shapefiles but they all differed significantly from the shape shown in the upper graph. Here some examples I've already tried but failed:

#Example 1
> gps_kde_raster<-raster(gps_kde)
> plot(gps_kde_raster)
> gps_kde_raster_95<-raster.vol(gps_kde_raster,p=0.95)
> plot(gps_kde_raster_95)
> writeRaster(gps_kde_raster_95,'gps.tif',options=c('TFW=YES'))
> gps_kde_95_polygon<-rasterToPolygons(gps_kde_raster_95,dissolve=TRUE)
> plot(gps_kde_95_polygon)
> writeOGR(gps_kde_95_polygon, ".", "gps_kde_95", driver="ESRI Shapefile")

Result: Raster layer with KDE-95%-area larger than the original (probably artefact of the low resolution)

#Example 2
> image(gps_kde$eval.points[[1]],gps_kde$eval.points[[2]],gps_kde$estimate)
> gps_kde_image<-image2Grid(list(x=gps_kde$eval.points[[1]],y=gps_kde$eval.points[[2]],z=gps_kde$estimate))
> contour(gps_kde_image,add=TRUE)

Result: Raster layer with KDE-area smaller than the original

#Example 3
> gps_raster<-raster(gps_kde)
> plot(gps_raster)
> gps_raster_max<-cellStats(gps_raster,max)
> gps_raster_final<-calc(gps_raster,fun=function(x){x/gps_raster_max})
> plot(gps_raster_final)
> writeRaster(gps_raster_final,'gps.tif',options=c('TFW=YES'))

Result: Raster layer with 0-1 probabilites as z-variable. Good, but low resolution leads to unreliable contour-creation

#Example 4
> r<-raster(extent(440000,465000,5655000,5680000),nrows=nrow(gps_kde$estimate),ncols=ncol(gps_kde$estimate))
> r[]<-gps_kde$estimate
> r<-flip(r,direction='y')
> rcont<-rasterToContour(r,levels=gps_kde$cont)
> (cont.values<-gps_kde$cont[grep(paste(c("50","75","95"),collapse="|"),names(gps_kde$cont))])
> rcont.gt50<-rcont[which(rcont@data[,1] %in% cont.values),]
> plot(rcont.gt50, main="50%, 75% and 95% volume contours")

Result: Vector layer with totally different KDE-boundaries than the original

#Example 5
> hts<-contourLevels(gps_kde, prob = c(0.5, 0.75, 0.95))
> c50<-contourLines(gps_kde$eval.points[[1]],gps_kde$eval.points[[2]],gps_kde$estimate,level=hts[1])
> c75<-contourLines(gps_kde$eval.points[[1]],gps_kde$eval.points[[2]],gps_kde$estimate,level=hts[2])
> c95<-contourLines(gps_kde$eval.points[[1]],gps_kde$eval.points[[2]],gps_kde$estimate,level=hts[3])
> Polyc50<-Polygon(c50[[1]][-1])
> Polysc50<-Polygons(list(Polyc50),"c50")
> Polyc75<-Polygon(c75[[1]][-1])
> Polysc75<-Polygons(list(Polyc75),"c75")
> Polyc95<-Polygon(c95[[1]][-1])
> Polysc95<-Polygons(list(Polyc95),"c95")
> spObj<-SpatialPolygons(list(Polysc50,Polysc75,Polysc95),1:3)
> axu.df<-data.frame(ID=1:length(spObj))
> rownames(axu.df)<-c("c50","c75","c95")
> spDFObj<-SpatialPolygonsDataFrame(spObj, axu.df)
> plot(spDFObj)

Result: Vector layer with totally different KDE-boundaries than the original

As you see, there many different approaches circulating, but none of them are scientifically useful. So is there anybody out there who knows how to export the real contour lines of my KDE (50%, 75% and 95%) as a polygon or polyline shapefile so I can load it into QGIS for further analysis (same boundaries as shown in second graph)?

Re: Save/Export specific contour levels for a KDE (package 'ks' in R) as ESRI shapefile

Verfasst: Di Sep 03, 2019 9:48 pm
von EDi
Bitte ein reproduzierbares Beispiel posten.

Re: Save/Export specific contour levels for a KDE (package 'ks' in R) as ESRI shapefile

Verfasst: Mi Sep 04, 2019 8:21 am
von raoulreding
Anbei ein reduzierter Originaldatensatz.

Re: Save/Export specific contour levels for a KDE (package 'ks' in R) as ESRI shapefile

Verfasst: Mi Okt 16, 2019 2:32 pm
von raoulreding
Das Problem konnte ich bislang leider noch nicht lösen. Konnte ziwschenzeitlich jemand anders schon mehr in Erfahrung bringen?