# blocked cv using ::blockCV library(blockCV) mat1 <- as.matrix(cbind(map_dat$LONG,map_dat$LAT,pred_dat$t_v)) mycrs <- "+proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs" rast1 <- raster::raster(mat1, crs = mycrs) spatialAutoRange(rasterLayer = rast1,showPlots = T) #crs(rast1) <- "+proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs" library(sp) # spatialBlock(speciesData = FirstPoints) # optionally optimize blocks # sac <- spatialAutoRange(rasterLayer = rast1, # sampleNumber = 5000, # doParallel = TRUE, # showPlots = TRUE) FirstPoints <- SpatialPoints(coords = cbind(map_dat$LAT,map_dat$LONG,pred_dat$t_v)) sb <- spatialBlock(speciesData = FirstPoints, theRange = 341126, # size of the blocks in meters k = 5, # number of folds selection = "random", iteration = 100, # find evenly dispersed folds biomod2Format = TRUE) # these are our blocked folds length(sb$foldID) sb$plots # cross validation with spatial blocks as folds k <- 5 D2 <- function(model) 1-model$deviance/model$null.deviance ell <- function(obs, pred) sum(log(dbinom(obs, prob=pred, size=1))) rmse <- function(actual, predicted) sqrt(mean((actual - predicted)^2)) rmse(100,1) pred_dat <- cbind(pred_dat, "foldid"=sb$foldID) mod1_loli <- 1:k mod2_loli <- 1:k mod1_auc <- 1:k mod2_auc <- 1:k mod1_rmse <- 1:k mod2_rmse <- 1:k predlist1 <- list() modlist1 <- list() predlist2 <- list() modlist2 <- list() # for auc library(Metrics) # cv of glm for (i in 1:k) { mod_dat <- pred_dat %>% filter(foldid!=i) val_dat <- pred_dat %>% filter(foldid==i) fmt <- glm(t_v~PRE_SUMMER+T_SUMMER, data=mod_dat, family = "binomial") # replace with fmt <- yourmodelobject$call preds1 <- plogis(predict.glm(fmt, newdata=val_dat)) mod1_rmse[i] <- rmse(actual = val_dat$t_v, predicted = preds1) mod1_loli[i] <- ell(obs = val_dat$t_v, pred = preds1) mod1_auc[i] <- auc(actual = val_dat$t_v, predicted = preds1) modlist1[[i]] <- fmt predlist1[[i]] <- preds1 } # cv of random forest library(randomForest) for (i in 1:k) { mod_dat <- pred_dat %>% filter(foldid!=i) val_dat <- pred_dat %>% filter(foldid==i) fmt <- randomForest(t_v~PRE_SUMMER+T_SUMMER, data=mod_dat, type="response") # replace with fmt <- yourmodelobject$call preds2 <- predict(fmt, newdata=val_dat) mod2_rmse[i] <- rmse(actual = val_dat$t_v, predicted = preds2) mod2_loli[i] <- ell(obs = val_dat$t_v, pred = preds2) # still can't get log lik with the rf mod2_auc[i] <- auc(actual = val_dat$t_v, predicted = preds2) modlist2[[i]] <- fmt predlist2[[i]] <- preds2 }