Commit a37bca88 authored by Glad Anouk's avatar Glad Anouk
Browse files

Revision

parent 98aea16e
No related merge requests found
Showing with 195 additions and 28 deletions
+195 -28
No preview for this file type
No preview for this file type
No preview for this file type
No preview for this file type
...@@ -27,6 +27,9 @@ Rpackage_path<-"/home/anouk/R/x86_64-pc-linux-gnu-library/3.2/" ...@@ -27,6 +27,9 @@ Rpackage_path<-"/home/anouk/R/x86_64-pc-linux-gnu-library/3.2/"
``` ```
``` {r echo=TRUE, eval=TRUE, tidy=FALSE, message=FALSE, results='hide'} ``` {r echo=TRUE, eval=TRUE, tidy=FALSE, message=FALSE, results='hide'}
Rpackage_path<-"/home/anouk/R/x86_64-pc-linux-gnu-library/3.2/"
.libPaths(Rpackage_path)
library(shapefiles) library(shapefiles)
library(logspline) library(logspline)
library(rgdal) library(rgdal)
...@@ -54,6 +57,8 @@ library(grid) ...@@ -54,6 +57,8 @@ library(grid)
library(gtable) library(gtable)
library(gridExtra) library(gridExtra)
library(RStoolbox) library(RStoolbox)
library(DSsim)
library(cvAUC)
##### linux paths ###### ##### linux paths ######
data_path<-"/media/anouk/DATAPART1/These/Data/DATAGTJPROPRE/virtual_sp/code_v_sp/Spatial_sampling_bias/Data/" data_path<-"/media/anouk/DATAPART1/These/Data/DATAGTJPROPRE/virtual_sp/code_v_sp/Spatial_sampling_bias/Data/"
...@@ -616,7 +621,7 @@ plot(region, plot.units = "m") ...@@ -616,7 +621,7 @@ plot(region, plot.units = "m")
design <- make.design(transect.type = "line", design <- make.design(transect.type = "line",
design.details = c("Parallel", "Systematic"), design.details = c("Parallel", "Systematic"),
region.obj = region, region.obj = region,
spacing = 300) spacing = 400)
transects <- generate.transects(design, region = region) transects <- generate.transects(design, region = region)
plot(region, plot.units = "m") plot(region, plot.units = "m")
...@@ -664,7 +669,7 @@ writeOGR(transect4, dsn = '.', layer = "Regular_transect", driver = "ESRI Shapef ...@@ -664,7 +669,7 @@ writeOGR(transect4, dsn = '.', layer = "Regular_transect", driver = "ESRI Shapef
transect<-shapefile(paste0(data_path, 'Regular_transect.shp')) transect<-shapefile(paste0(data_path, 'Regular_transect.shp'))
#Use raster of environmental variable in order set the right extent #Use raster of environmental variable in order set the right extent
a<-raster(extent(new_raster3)) #choose a raster to set up extent a<-raster(extent(new_env_raster)) #choose a raster to set up extent
res(a)<-c(25,25) res(a)<-c(25,25)
#Create a raster representing trails (same extent and resolution as environmental variables rasters) #Create a raster representing trails (same extent and resolution as environmental variables rasters)
transectR <- rasterize(transect, a) transectR <- rasterize(transect, a)
...@@ -4568,19 +4573,26 @@ row_spec(4:6, color = "black", background = "#fff") ...@@ -4568,19 +4573,26 @@ row_spec(4:6, color = "black", background = "#fff")
*Open environmental data* *Open environmental data*
```{r echo=TRUE, eval=TRUE, tidy=TRUE, message=FALSE} ```{r echo=TRUE, eval=TRUE, tidy=TRUE, message=FALSE}
new_env_raster <- stack( new_env_raster <- stack(paste0(data_path, "JURAraster_MV_1_8ha_Mean_H.nb10_20relative_density.tif"),
paste0(data_path, "JURAraster_H.nb10_20relative_density.tif"),
paste0(data_path, "JURAraster_H.nb20_30relative_density.tif"), paste0(data_path, "JURAraster_MV_1_8ha_Mean_H.nb2_5ratio.tif"),
paste0(data_path, "JURAraster_H.simpson.tif"),
paste0(data_path, "JURAraster_MV_1_8ha_Mean_H.simpson.tif"),
paste0(data_path, "JURAraster_MV_1_8ha_Mean_H.p25.tif"),
paste0(data_path, "JURAraster_MV_1_8ha_sd_H.nb20_30relative_density.tif"),
paste0(data_path, "JURAraster_MV_1_8ha_sd_H.nb2_5ratio.tif"),
quick=TRUE quick=TRUE
) )
names(new_env_raster) <- c("Canopy1020", "Canopy2030", "Simpson") names(new_env_raster) <- c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
#crs(new_env_raster)<-"+proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +vunits=m +no_defs" crs(new_env_raster)<-"+proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +vunits=m +no_defs"
new_env_raster<-normImage(new_env_raster) new_env_raster<-normImage(new_env_raster)
...@@ -4592,7 +4604,7 @@ new_env_raster_tot <- overlay(new_env_raster, rpoly, fun = function(x, y) { ...@@ -4592,7 +4604,7 @@ new_env_raster_tot <- overlay(new_env_raster, rpoly, fun = function(x, y) {
x[is.na(y[])] <- NA x[is.na(y[])] <- NA
return(x) return(x)
}) })
names(new_env_raster_tot) <-c( "Canopy1020", "Canopy2030", "Simpson" ) names(new_env_raster_tot) <-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
predictor_tot<-rasterToPoints(new_env_raster_tot) predictor_tot<-rasterToPoints(new_env_raster_tot)
predictor_tot<-na.omit(predictor_tot) predictor_tot<-na.omit(predictor_tot)
...@@ -4650,7 +4662,7 @@ new_env_raster_S <- overlay(new_env_raster, rpoly_S, fun = function(x, y) { ...@@ -4650,7 +4662,7 @@ new_env_raster_S <- overlay(new_env_raster, rpoly_S, fun = function(x, y) {
x[is.na(y[])] <- NA x[is.na(y[])] <- NA
return(x) return(x)
}) })
names(new_env_raster_S) <-c( "Canopy1020", "Canopy2030", "Simpson" ) names(new_env_raster_S) <-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
# Non Sytematic # Non Sytematic
...@@ -4663,7 +4675,7 @@ new_env_raster_R <- overlay(new_env_raster, rpoly_R, fun = function(x, y) { ...@@ -4663,7 +4675,7 @@ new_env_raster_R <- overlay(new_env_raster, rpoly_R, fun = function(x, y) {
return(x) return(x)
}) })
names(new_env_raster_R) <-c( "Canopy1020", "Canopy2030", "Simpson") names(new_env_raster_R) <-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
``` ```
*Select data (systemtic/Conventional design and Female/male)* *Select data (systemtic/Conventional design and Female/male)*
```{r echo=TRUE, eval=TRUE, tidy=TRUE, message=FALSE} ```{r echo=TRUE, eval=TRUE, tidy=TRUE, message=FALSE}
...@@ -4745,6 +4757,92 @@ pres_pointsGDT<- data.frame(x=OBSGDTshp$coords.x1, y=OBSGDTshp$coords.x2, presen ...@@ -4745,6 +4757,92 @@ pres_pointsGDT<- data.frame(x=OBSGDTshp$coords.x1, y=OBSGDTshp$coords.x2, presen
``` ```
```{r echo=TRUE, eval=TRUE, tidy=TRUE, message=FALSE}
#fonction for fold and AUC
iid <- function(data, V, predictors_vector){
.cvFolds <- function(Y, V){ #Create CV folds (stratify by outcome)
Y0 <- split(sample(which(Y==0)), rep(1:V, length=length(which(Y==0))))
Y1 <- split(sample(which(Y==1)), rep(1:V, length=length(which(Y==1))))
folds <- vector("list", length=V)
for (v in seq(V)) {folds[[v]] <- c(Y0[[v]], Y1[[v]])}
return(folds)
}
.doFit <- function(v, folds, data){ #Train/test glm for each fold
fit <- maxnet(as.vector(data[-folds[[v]],"presence"]), data[-folds[[v]], predictors_vector], f = FORMULA, regmult = 1)
pred <- predict(fit, newdata=data[folds[[v]], predictors_vector], type=c("exponential"), clamp=FALSE)
return(pred)
}
set.seed(158)
folds <- .cvFolds(Y=data$presence, V=V) #Create folds
predictions <- unlist(sapply(seq(V), .doFit, folds=folds, data=data)) #CV train/predict
predictions[unlist(folds)] <- predictions #Re-order pred va
# Get CV AUC and confidence interval
out <- ci.cvAUC(predictions=predictions, labels=data$presence, folds=folds, confidence=0.95)
#Get models for output
# list_fit<-list()
#
# list_fit <- lapply(folds, function(x) maxnet(as.vector(data[-x,"presence"]), data[-x, predictors_vector], f = FORMULA, regmult = 1) )
#
# OUT<-list(out, list_fit)
# return(OUT)
return(out)
}
iid2 <- function(data, V, predictors_vector){
.cvFolds <- function(Y, V){ #Create CV folds (stratify by outcome)
Y0 <- split(sample(which(Y==0)), rep(1:V, length=length(which(Y==0))))
Y1 <- split(sample(which(Y==1)), rep(1:V, length=length(which(Y==1))))
folds <- vector("list", length=V)
for (v in seq(V)) {folds[[v]] <- c(Y0[[v]], Y1[[v]])}
return(folds)
}
.doFit <- function(v, folds, data){ #Train/test glm for each fold
fit <- maxnet(as.vector(data[-folds[[v]],"presence"]), data[-folds[[v]], predictors_vector], f = FORMULA, regmult = 1)
data[, 10][data[, 10] != 0] <- 0
pred <- predict(fit, newdata=data[folds[[v]], predictors_vector], type=c("exponential"), clamp=FALSE)
return(pred)
}
set.seed(158)
folds <- .cvFolds(Y=data$presence, V=V) #Create folds
predictions <- unlist(sapply(seq(V), .doFit, folds=folds, data=data)) #CV train/predict
predictions[unlist(folds)] <- predictions #Re-order pred va
# Get CV AUC and confidence interval
out <- ci.cvAUC(predictions=predictions, labels=data$presence, folds=folds, confidence=0.95)
#Get models for output
# list_fit<-list()
#
# list_fit <- lapply(folds, function(x) maxnet(as.vector(data[-x,"presence"]), data[-x, predictors_vector], f = FORMULA, regmult = 1) )
#
# OUT<-list(out, list_fit)
# return(OUT)
return(out)
}
```
```{r echo=TRUE, eval=TRUE, tidy=TRUE, message=FALSE} ```{r echo=TRUE, eval=TRUE, tidy=TRUE, message=FALSE}
...@@ -4773,11 +4871,16 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -4773,11 +4871,16 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
out_F_1 <- iid(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
predictor<-rasterToPoints(new_env_raster_S) predictor<-rasterToPoints(new_env_raster_S)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
...@@ -4832,11 +4935,16 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -4832,11 +4935,16 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
out_F_2 <- iid(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
predictor<-rasterToPoints(new_env_raster_S) predictor<-rasterToPoints(new_env_raster_S)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
...@@ -4924,10 +5032,16 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -4924,10 +5032,16 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson", "SystTracks_distance5_5")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd", "SystTracks_distance5_5")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd", "SystTracks_distance5_5")
out_F_3 <- iid2(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
#Set all values to 0 #Set all values to 0
raster_tracks_syst_1<-raster_tracks_syst raster_tracks_syst_1<-raster_tracks_syst
...@@ -4943,7 +5057,7 @@ predictor$predicted<-predict(model_maxnet, predictor,type=c("exponential"), clam ...@@ -4943,7 +5057,7 @@ predictor$predicted<-predict(model_maxnet, predictor,type=c("exponential"), clam
maps_predicted_F_3<-rasterFromXYZ(predictor[, c("x","y","predicted")], res=c(25,25)) maps_predicted_F_3<-rasterFromXYZ(predictor[, c("x","y","predicted")], res=c(25,25))
colnames(predictor_tot2) <- c("x", "y", "SystTracks_distance5_5", "Canopy1020", "Canopy2030", "Simpson") colnames(predictor_tot2) <- c("x", "y", "SystTracks_distance5_5", "Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
predictor_tot2$predicted_tot2<-predict(model_maxnet, predictor_tot2,type=c("exponential"), clamp=FALSE) predictor_tot2$predicted_tot2<-predict(model_maxnet, predictor_tot2,type=c("exponential"), clamp=FALSE)
maps_predicted_F_3b<-rasterFromXYZ(predictor_tot2[, c("x","y","predicted_tot2")], res=c(25,25)) maps_predicted_F_3b<-rasterFromXYZ(predictor_tot2[, c("x","y","predicted_tot2")], res=c(25,25))
...@@ -4993,10 +5107,16 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -4993,10 +5107,16 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
out_F_4 <- iid(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
predictor<-rasterToPoints(new_env_raster_R) predictor<-rasterToPoints(new_env_raster_R)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
...@@ -5047,11 +5167,16 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -5047,11 +5167,16 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
out_F_5 <- iid(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
predictor<-rasterToPoints(new_env_raster_R) predictor<-rasterToPoints(new_env_raster_R)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
...@@ -5135,10 +5260,17 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -5135,10 +5260,17 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson", "BiaisTracks_distance5_5")]
predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd", "BiaisTracks_distance5_5")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd", "BiaisTracks_distance5_5")
out_F_6 <- iid2(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
#Set all values to 1 #Set all values to 1
raster_tracks_bias_1<-raster_tracks_bias raster_tracks_bias_1<-raster_tracks_bias
...@@ -5152,7 +5284,7 @@ predictor<-as.data.frame(predictor) ...@@ -5152,7 +5284,7 @@ predictor<-as.data.frame(predictor)
predictor$predicted<-predict(model_maxnet, predictor,type=c("exponential"), clamp=FALSE) predictor$predicted<-predict(model_maxnet, predictor,type=c("exponential"), clamp=FALSE)
colnames(predictor_tot2) <- c("x", "y", "BiaisTracks_distance5_5", "Canopy1020", "Canopy2030", "Simpson") colnames(predictor_tot2) <- c("x", "y", "BiaisTracks_distance5_5", "Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
predictor_tot2$predicted_tot2<-predict(model_maxnet, predictor_tot2,type=c("exponential"), clamp=FALSE) predictor_tot2$predicted_tot2<-predict(model_maxnet, predictor_tot2,type=c("exponential"), clamp=FALSE)
...@@ -5211,10 +5343,16 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -5211,10 +5343,16 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
out_M_1 <- iid(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
predictor<-rasterToPoints(new_env_raster_S) predictor<-rasterToPoints(new_env_raster_S)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
...@@ -5261,11 +5399,16 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -5261,11 +5399,16 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
out_M_2 <- iid(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
predictor<-rasterToPoints(new_env_raster_S) predictor<-rasterToPoints(new_env_raster_S)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
...@@ -5312,10 +5455,18 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -5312,10 +5455,18 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson", "SystTracks_distance5_5")] model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd", "SystTracks_distance5_5")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd", "SystTracks_distance5_5")
out_F_3 <- iid2(data=pres_bgGDT, V=10, var_names)
predictor<-rasterToPoints(new_env_rastersyst3) predictor<-rasterToPoints(new_env_rastersyst3)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
predictor<-as.data.frame(predictor) predictor<-as.data.frame(predictor)
...@@ -5323,7 +5474,7 @@ predictor<-as.data.frame(predictor) ...@@ -5323,7 +5474,7 @@ predictor<-as.data.frame(predictor)
predictor$predicted<-predict(model_maxnet, predictor,type=c("exponential"), clamp=FALSE) predictor$predicted<-predict(model_maxnet, predictor,type=c("exponential"), clamp=FALSE)
#predictor_tot2<-predictor_tot2[, c(-7,-8)] #predictor_tot2<-predictor_tot2[, c(-7,-8)]
colnames(predictor_tot2) <- c("x", "y", "SystTracks_distance5_5", "Canopy1020", "Canopy2030", "Simpson") colnames(predictor_tot2) <- c("x", "y", "SystTracks_distance5_5", "Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
predictor_tot2$predicted_tot2<-predict(model_maxnet, predictor_tot2,type=c("exponential"), clamp=FALSE) predictor_tot2$predicted_tot2<-predict(model_maxnet, predictor_tot2,type=c("exponential"), clamp=FALSE)
...@@ -5375,11 +5526,16 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -5375,11 +5526,16 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
out_M_4 <- iid(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
predictor<-rasterToPoints(new_env_raster_R) predictor<-rasterToPoints(new_env_raster_R)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
...@@ -5427,11 +5583,16 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -5427,11 +5583,16 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
out_M_5 <- iid(data=pres_bgGDT, V=10, var_names)
#use it to predict results #use it to predict results
predictor<-rasterToPoints(new_env_raster_R) predictor<-rasterToPoints(new_env_raster_R)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
...@@ -5476,17 +5637,23 @@ pres_bgGDT <- na.omit(pres_bgGDT) ...@@ -5476,17 +5637,23 @@ pres_bgGDT <- na.omit(pres_bgGDT)
Pres_data<-as.vector(pres_bgGDT$presence) Pres_data<-as.vector(pres_bgGDT$presence)
predictor_maxnet<-pres_bgGDT[, c( "Canopy1020", "Canopy2030", "Simpson", "BiaisTracks_distance5_5")] predictor_maxnet<-pres_bgGDT[, c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd", "BiaisTracks_distance5_5")]
FORMULA<-maxnet.formula(Pres_data, predictor_maxnet, classes = "q")
model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1) model_maxnet<-maxnet(Pres_data, predictor_maxnet, f = maxnet.formula(Pres_data, predictor_maxnet, classes = "q"), regmult = 1)
var_names<-c("Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd", "BiaisTracks_distance5_5")
out_F_6 <- iid2(data=pres_bgGDT, V=10, var_names)
predictor<-rasterToPoints(new_env_rasterbias3) predictor<-rasterToPoints(new_env_rasterbias3)
predictor<-na.omit(predictor) predictor<-na.omit(predictor)
predictor<-as.data.frame(predictor) predictor<-as.data.frame(predictor)
predictor$predicted<-predict(model_maxnet, predictor,type=c("exponential"), clamp=FALSE) predictor$predicted<-predict(model_maxnet, predictor,type=c("exponential"), clamp=FALSE)
colnames(predictor_tot2) <- c("x", "y", "BiaisTracks_distance5_5", "Canopy1020", "Canopy2030", "Simpson") colnames(predictor_tot2) <- c("x", "y", "BiaisTracks_distance5_5", "Canopy1020","penetrationratio0205","Simpson", "Q25", "Canopy2030sd","penetrationratio0205sd")
predictor_tot2$predicted_tot2<-predict(model_maxnet, predictor_tot2,type=c("exponential"), clamp=FALSE) predictor_tot2$predicted_tot2<-predict(model_maxnet, predictor_tot2,type=c("exponential"), clamp=FALSE)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment