0

所以我一直在尝试制作一个未来的 Salsola Tragus 物种的投影模型,但我在接近尾声时遇到了一个障碍

 bc <- bioclim(presvals[,c('wc2.1_10m_bioc_CNRM-CM6-1_ssp370_2021-2040','-wc2.1_10m_prec_CNRM-CM6-1_ssp370_2021-2040','wc2.1_10m_tmax_CNRM-CM6-1_ssp370_2021-2040','wc2.1_10m_tmin_CNRM-CM6-1_ssp370_2021-2040')])
Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 'bioclim': subscript out of bounds 

我一直在搞乱配对功能,但我无法弄清楚其他任何事情。

下面是我一直在运行的代码

setwd("~/gfds/Salsola Files")
#loading in the needed library
library(dismo)
library(maptools)
library(raster)
library(sp)
library(rgeos)
library(caret)
library(rJava)
library(gdalUtils)
library(rgdal)

#---misc variables manually defined-----------------------------------------------------------------

e <- extent (-180,180, -80, 80) #Set extent of area (long_min, long_max, lat_min, lat_max)
n_abs <- 1000 

datafile <- "Salsola Tragus Data.csv"

data <- read.csv(datafile,sep = ",", quote = "\"",dec = ".", fill = TRUE, comment.char = "")

data <- data[,c("Scientific.Names","Longitude","Latitude")]

salsola.coords <- data[,c("Longitude", "Latitude")]

data.sdf<-data              #spacial data frame for salsola data

coordinates(data.sdf) <- ~Longitude+Latitude

projection(data.sdf) <- CRS('+proj=longlat +datum=WGS84')

#---Load worldclim data
#---Read in the Tif Files -------------------------------------------------------------
#Load bioclim stack

future1<-'wc2.1_10m_bioc_CNRM-CM6-1_ssp370_2021-2040.tif' 
future2<-'wc2.1_10m_prec_CNRM-CM6-1_ssp370_2021-2040.tif'
future3<-'wc2.1_10m_tmax_CNRM-CM6-1_ssp370_2021-2040.tif'
future4<-'wc2.1_10m_tmin_CNRM-CM6-1_ssp370_2021-2040.tif' 

#----------------creating the rasters and stacking them on top of eachother ----------------
ras1<-raster(future1)
ras2<-raster(future2)
ras3<-raster(future3)
ras4<-raster (future4)

Rastack<-stack(ras1,ras2,ras3,ras4) #stacking up all the rasters so that htey can replaced the layered one 

worldclim.extent <-stack(crop(Rastack,e)) # adding hte extent to the stack so that hte raster is trimmed to the right size 

#---pseudo absence simulation------------------------------------------------------------------------------------------------------------------
#Define a mask from worldclim data
mask_land <- raster(worldclim.extent@layers[[1]])

# Create circles with a radius of 100 km
x <- circles(data.sdf, d=10000, lonlat=TRUE)
presvals <- extract(worldclim.extent, salsola.coords) # start here --------------------------------

# setting random seed to always create the same
# random set of points for this example

set.seed(0)

stvals <- extract(worldclim.extent, randomPoints(worldclim.extent, n_abs)) #simulated absence values

pb <- c(rep(1, nrow(presvals)), rep(0, nrow(stvals)))#presence & absence binary array

sdmdata <- data.frame(cbind(pb, rbind(presvals, stvals)))#dataframe of presence and absence data

#Plot using "pairs" function--------------------------------------------------------------------------------------------------------------------------------
pairs(sdmdata[,2:5], cex=0.1)

#Save outputs to be used in next section
saveRDS(sdmdata, "s_tragus.Rds")
saveRDS(presvals, "s_tragus_pvals.Rds")

#---Model Fitting
#GLM: Linear Models

m2 = glm(pb ~ ., data=sdmdata)

#MaxEnt model
#get maxent() #Run this prior to attempting to load dismo

#Create training data for maxent------------------------------------------------------------------------------------------------------------------------------
group<-kfold(salsola.coords,5)

pres_test <- salsola.coords[group==1, ]

pres_train <-salsola.coords[group !=1,]

#Run MaxEnt: I am not totally sure how "factors" works, ------------------------------------------------------------------------------------------------------
#but it tells MaxEnt which layers we want to consider "categorical". I've set
#this series of model creation steps up to try each of our chosen variables, 
#as well as one which doesn't specify a factor (max).

#FYI, MaxEnt has a "removeDuplicates" Boolean option which will remove points
#which fall in the same grid cell. That might be really useful if your abiotic
#rasters have varying resolution.
maxcom <- maxent(worldclim.extent,pres_train)#####THIS SHOULD BE PRES_TRAIN####

maxcom2 <- maxent(worldclim.extent,pres_train,factors='wc2.1_10m_bioc_CNRM-CM6-1_ssp370_2021-2040')
maxcom3 <- maxent(worldclim.extent,pres_train,factors='wc2.1_10m_prec_CNRM-CM6-1_ssp370_2021-2040')
maxcom4 <- maxent(worldclim.extent,pres_train,factors='wc2.1_10m_tmax_CNRM-CM6-1_ssp370_2021-2040')
maxcom5 <- maxent(worldclim.extent,pres_train,factors='wc2.1_10m_tmin_CNRM-CM6-1_ssp370_2021-2040')

#bioclim model-------------------------------------------------------------------------------------------------------------------------------------------
bc <- bioclim(presvals[,c('wc2.1_10m_bioc_CNRM-CM6-1_ssp370_2021-2040','-wc2.1_10m_prec_CNRM-CM6-1_ssp370_2021-2040','wc2.1_10m_tmax_CNRM-CM6-1_ssp370_2021-2040','wc2.1_10m_tmin_CNRM-CM6-1_ssp370_2021-2040')])

任何能指引我正确方向的东西都会有很大帮助

4

0 回答 0