-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathfor loop dismo rcp45 both.R
51 lines (46 loc) · 2.55 KB
/
for loop dismo rcp45 both.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
library(raster)
library(rgdal)
library(dismo)
library(sp)
library(rJava)
library(readxl)
library(mapview)
library(ggplot2)
library(jsonlite)
require(sf)
library(devtools)
library(xlsx)
library(dplyr)
# options(java.parameters = "-Xmx8000m")
#latlong <- CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0")
#SR.ORG8287 <- CRS('+proj=cea +lon_0=0 +lat_ts=0 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs')
#extent.asia.sr <- extent(11000000, 13000000, 900000, 2600000)
# Current
list.raster.full <- list.files("F:/Working/2018/PhD_research/enviromental_variables/current", full.names = T, pattern = ".tif")
predictors <- stack(list.raster.full)
predictors <- subset(predictors, c("bio2","bio10","bio11","bio12","bio18","bio19","karst","forest","grassland","farmland","urban")) #select variables
# 2050
#rcp45
land.cover.2050.rcp45 <- list.files("F:/Working/2018/PhD_research/enviromental_variables/2050/rcp45/land-cover", full.names = T, pattern = ".tif")
land.cover.2050.rcp45 <- stack(land.cover.2050.rcp45)
karst <- raster('F:/Working/2018/PhD_research/enviromental_variables/current/karst.tif')
#import species list
species.list <- read_excel('F:/Working/2018/PhD_research/SDM output R/17_5_23_BC/list of species.xlsx')
species.list <- species.list%>%filter(Num_occurence_all > 19)
nrow(species.list)
#predict future under combined change
for(sp in species.list$Scientific_name){
load(paste0('F:/Working/2018/PhD_research/SDM output R/17_5_23_BC/response/',sp,'.RData'))
maxent.result <- read.csv(paste0('F:/Working/2018/PhD_research/SDM output R/17_5_23_BC/response/',sp,'/maxentResults.csv'), header=T, sep=',')
for(gcm in c('ac', 'bc','ca','cm','cn','cs','gf','gi','ip','mi')){
climate.2050.rcp45 <- list.files(paste0("F:/Working/2018/PhD_research/enviromental_variables/2050/rcp45/",gcm), full.names = T, pattern = ".tif")
climate.2050.rcp45 <- stack(climate.2050.rcp45)
climate.2050.rcp45 <- subset(climate.2050.rcp45, c("bio2","bio10","bio11","bio12","bio18","bio19"))
predictors.2050.rcp45 <- stack(climate.2050.rcp45, land.cover.2050.rcp45, karst)
px.2050.rcp45.both <- predict(predictors.2050.rcp45, xm, ext=extent(11000000, 13000000, 900000, 2600000))
species.2050.rcp45.both.bin <- px.2050.rcp45.both > maxent.result$X10.percentile.training.presence.Cloglog.threshold
writeRaster(species.2050.rcp45.both.bin, filename = paste0('F:/Working/2018/PhD_research/SDM output R/17_5_23_BC/output_maps/',sp,' 2050 rcp45 both ',gcm,'.tif'),format='GTiff', overwrite=T)
print(gcm)
}
print(sp)
}