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

... zu anderer statistischer Software, zu Datenbanken und Programmiersprachen.

Moderatoren: EDi, jogo

Antworten
raoulreding
Beiträge: 3
Registriert: Di Sep 03, 2019 1:55 pm

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

Beitrag 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) 2612 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) 2612 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)?
Benutzeravatar
EDi
Beiträge: 1599
Registriert: Sa Okt 08, 2016 3:39 pm

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

Beitrag von EDi »

Bitte ein reproduzierbares Beispiel posten.
Bitte immer ein reproduzierbares Minimalbeispiel angeben. Meinungen gehören mir und geben nicht die meines Brötchengebers wieder.

Dieser Beitrag ist lizensiert unter einer CC BY 4.0 Lizenz
Bild.
raoulreding
Beiträge: 3
Registriert: Di Sep 03, 2019 1:55 pm

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

Beitrag von raoulreding »

Anbei ein reduzierter Originaldatensatz.
Dateianhänge
gps.csv
(7.21 KiB) 143-mal heruntergeladen
raoulreding
Beiträge: 3
Registriert: Di Sep 03, 2019 1:55 pm

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

Beitrag von raoulreding »

Das Problem konnte ich bislang leider noch nicht lösen. Konnte ziwschenzeitlich jemand anders schon mehr in Erfahrung bringen?
Antworten