1

我想使用 BIEN 使用的相同方法为物种出现绘制多边形，所以我可以同时使用我的多边形和他们的多边形。他们使用 Maxent 来模拟物种的出现，当它们有更多的出现点时。

例如，这是一个 BIEN 多边形：

library(BIEN)
  Mormolyca_ringens<- BIEN_ranges_load_species(species = "Mormolyca ringens")
#And this is a polygon, yes. A SpatialPolygonsDataFrame. 
  plot(wrld_simpl, xlim=c(-100,-40), ylim=c(-30,30), axes=TRUE,col="light yellow", bg="light blue")
  plot(Mormolyca_ringens, col="green", add=TRUE)

Mormolyca ringens 多边形 Mormolyca ringens 多边形

好的，然后我正在尝试绘制我的多边形，因为 BIEN 缺少一些我需要的物种。

# first, you need to download the Maxent software here: http://biodiversityinformatics.amnh.org/open_source/maxent/
#and paste the "maxent.jar" file in the ’java’ folder of the ’dismo’ package, which is here:

system.file("java", package="dismo")

#You have to do this **before** loading the libraries

#install.packages("rJava")
library(rJava)
#If you get the message that cannot load this library, it's possible that your version of java is not 64bit. 
#Go to Oracle and install Java for windows 64bit.
#If library still doesn't load: Look in your computer for the path where the java's jre file is and paste in the code below
Sys.setenv(JAVA_HOME="your\\path\\for\\jre") #mine is "C:\\Program Files\\Java\\jre1.8.0_144", for example

library(rJava)
library(dismo)
library(maptools)

#Giving credits: I wrote the following code based on this tutorial: https://cran.r-project.org/web/packages/dismo/vignettes/sdm.pdf


#Preparing the example data - the map
data(wrld_simpl) 
ext = extent(-90, -32, -33, 23)

#Preparing the example data -  presence data for Bradypus variegatus
file <- paste(system.file(package="dismo"), "/ex/bradypus.csv", sep="")
bradypus <- read.table(file, header=TRUE, sep=',')
bradypus <- bradypus[,-1] #don't need th first col

#Getting the predictors (the variables)
files <- list.files(path=paste(system.file(package="dismo"),
                               '/ex', sep=''), pattern='grd', full.names=TRUE )
predictors <- stack(files)



#making a training and a testing set.
group <- kfold(bradypus, 5)
pres_train <- bradypus[group != 1, ]
pres_test <- bradypus[group == 1, ]

#Creating the background
backg <- randomPoints(predictors, n=1000, ext=ext, extf = 1.25)
colnames(backg) = c('lon', 'lat')
group <- kfold(backg, 5)
backg_train <- backg[group != 1, ]
backg_test <- backg[group == 1, ]


# Running maxent
  xm <- maxent(predictors, pres_train, factors='biome')
  plot(xm)
  

#A response plot:
  response(xm)

# Evaluating and predicting 
  e <- evaluate(pres_test, backg_test, xm, predictors)
   px <- predict(predictors, xm, ext=ext,  progress='text', overwrite=TRUE)
                
  #Checking result of the prediction
   par(mfrow=c(1,2))
  plot(px, main='Maxent, raw values')
  plot(wrld_simpl, add=TRUE, border='dark grey')
  tr <- threshold(e, 'spec_sens')
  plot(px > tr, main='presence/absence')
  plot(wrld_simpl, add=TRUE, border='dark grey')
  points(pres_train, pch='+')

此时，我有以下图像：

预测例子的发生 预测例子的发生

我正在尝试使用此代码从该栅格制作一个多边形：

predic_pol<-rasterToPolygons(px )

并且：

px_rec<-reclassify(px, rcl=0.5,  include.lowest=FALSE)
px_pol<-rasterToPolygons(px_rec)

但我不断得到我范围的像素版本

你能给我一个提示，这样我就可以从这个栅格中提取一个多边形，就像 BIEN 的一样？（我也是建模和 R 的新手......欢迎任何提示）

编辑：这是 px 控制台输出：

> px
class       : RasterLayer 
dimensions  : 172, 176, 30272  (nrow, ncol, ncell)
resolution  : 0.5, 0.5  (x, y)
extent      : -120, -32, -56, 30  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
data source : C:\Users\thai\Documents\ORCHIDACEAE\Ecologicos\w2\predictions\Trigonidiumobtusum_prediction.grd 
names       : layer 
values      : 6.705387e-06, 0.9999983  (min, max)

先感谢您

编辑2：解决方案

感谢@Val，我得到了这个：

#Getting only the values>tr to make the polygon 
#"tr" is what gives me the green raster instear of the multicolour one 
pol <- rasterToPolygons(px>tr,function(x) x == 1,dissolve=T)
#Ploting
plot(wrld_simpl, xlim=c(-120,-20), ylim=c(-60,10), axes=TRUE,col="light yellow", bg="light blue")
plot(pol, add=T, col="green")

现在我得到了我想要的！谢谢！（图中的多边形不一样，只是因为我在得到@Val 的答案时使用了我在环境中拥有的不同数据集）

奖金问题：

你知道如何平滑边缘，所以我得到一个非像素化的多边形吗？

4

1 回答 1

1

我不知道BIEN，所以我并没有真正看你示例的这一部分。我只是将您的问题/问题概括为以下内容：

您有一个二进制栅格（0 表示不存在，1 表示存在）并且您希望将所有具有 1 的区域转换为多边形。

至于您的栅格，您的值不是 0 和 1，而是更基本上是0 和基本上px是 1 ，这有点奇怪。但如果这是一个问题，那可以很容易解决。

因此，我尝试仅使用巴西地区重新创建您的示例：

library(raster)
library(rgeos)

# get Brasil borders

shp <- getData(country = 'BRA',level=0)

#create binary raster


r <- raster(extent(shp),resolution=c(0.5,0.5))
r[] <- NA # values have to be NA for the buffering

# take centroid of Brasil as center of species presence
cent <- gCentroid(shp)

# set to 1 
r[cellFromXY(r,cent)] <- 1

# buffer presence
r <- buffer(r,width=1000000)

# set rest 0
r[is.na(r)] <- 0

# mask by borders
r <- mask(r,shp)

我猜这与您的栅格足够接近：

存在图

所以现在转换为多边形：

pol <- rasterToPolygons(r,function(x) x == 1,dissolve=T)

我使用一个函数仅获取值为 1 的像素。此外，我将多边形分解为没有单个像素多边形，而是一个区域。请参阅rasterToPolygons其他选项。

现在将边界和新多边形绘制在一起：

plot(shp)
plot(pol,col='red',add=T)

多边形

这就是分布的多边形。这是控制台输出：

> pol
class       : SpatialPolygonsDataFrame 
features    : 1 
extent      : -62.98971, -43.48971, -20.23512, -1.735122  (xmin, xmax, ymin, ymax)
coord. ref. : NA 
variables   : 1
names       : layer 
min values  :     1 
max values  :     1

希望有帮助！

编辑：奖励答案

您必须清楚，多边形的像素化边界代表了数据的准确表示。因此，对此的任何更改都意味着精度的损失。现在，根据您的目的，这可能无关紧要。

有多种方法可以实现它，无论是在光栅侧进行分解和平滑/过滤等，还是在多边形侧，您可以将特定的过滤器应用于这样的多边形。

如果它纯粹是美学，您可以gSimplify从rgeos包装中尝试：

# adjust tol for smoothness
pol_sm <- gSimplify(pol,tol=0.5)

plot(pol)
lines(pol_sm,col='red',lwd=2)

在此处输入图像描述

于 2017-08-03T11:35:51.383 回答