Merge branch 'HEAD'

DESCRIPTION

@@ -1,8 +1,8 @@
 Package: ENMTools
 Type: Package
 Title: Analysis of Niche Evolution using Niche and Distribution Models
-Version: 1.0
-Date: 2020-07-22
+Version: 1.0.1
+Date: 2020-08-25
 Encoding: UTF-8
 Authors@R: c(person(given = "Dan", family = "Warren", role = c("aut", "cre"), email = "dan.l.warren@gmail.com", comment = c(ORCID = "0000-0002-8747-2451")),
             person(given = "Nick", family = "Matzke", role = "ctb", email = "nick.matzke@anu.edu.au"),
@@ -29,7 +29,8 @@ Imports:
 	sp,
 	spatstat,
 	ggpubr,
-	rgdal
+	rgdal,
+	forcats
 Depends:
 	raster,
 	dismo
@@ -47,4 +48,9 @@ Suggests:
 	ResourceSelection,
 	ppmlasso,
 	fields,
-	reshape2
+	reshape2,
+	vip,
+	pdp,
+	fastshap,
+	viridis,
+	progress

NAMESPACE

@@ -8,6 +8,7 @@ export(calc.B2)
 export(check.bg)
 export(check.clade)
 export(check.species)
+export(clamp.env)
 export(drop.species)
 export(enmtools.aoc)
 export(enmtools.bc)
@@ -24,6 +25,7 @@ export(enmtools.ppmlasso)
 export(enmtools.rf.ranger)
 export(enmtools.rf)
 export(enmtools.species)
+export(enmtools.vip)
 export(env.breadth)
 export(env.evaluate)
 export(env.overlap)
@@ -38,6 +40,7 @@ export(node.overlap)
 export(install.extras)
 export(interactive.plot.enmtools.species)
 export(interactive.plot.enmtools.model)
+export(interactive.plot)
 export(point.overlap)
 export(raster.cor.plot)
 export(geog.range.overlap)
@@ -130,47 +133,15 @@ importFrom(dismo,evaluate)
 importFrom(dismo,bioclim)
 importFrom(dismo,maxent)
 importFrom(ENMeval,get.block)
+importFrom(forcats,fct_reorder)
 importFrom(graphics,plot.new)
 importFrom(graphics,barplot)
 importFrom(graphics,abline)
 importFrom(graphics,plot)
 importFrom(graphics,points)
 importFrom(graphics,title)
 importFrom(graphics,par)
-importFrom(ggplot2,geom_abline)
-importFrom(ggplot2,geom_contour)
-importFrom(ggplot2,geom_point)
-importFrom(ggplot2,geom_raster)
-importFrom(ggplot2,geom_vline)
-importFrom(ggplot2,ggplot)
-importFrom(ggplot2,ggtitle)
-importFrom(ggplot2,labs)
-importFrom(ggplot2,layer)
-importFrom(ggplot2,qplot)
-importFrom(ggplot2,theme_bw)
-importFrom(ggplot2,theme_classic)
-importFrom(ggplot2,xlab)
-importFrom(ggplot2,xlim)
-importFrom(ggplot2,ylim)
-importFrom(ggplot2,geom_tile)
-importFrom(ggplot2,coord_fixed)
-importFrom(ggplot2,theme)
-importFrom(ggplot2,element_text)
-importFrom(ggplot2,theme)
-importFrom(ggplot2,stat_density_2d)
-importFrom(ggplot2,qplot)
-importFrom(ggplot2,guides)
-importFrom(ggplot2,guide_colourbar)
-importFrom(ggplot2,aes)
-importFrom(ggplot2,aes_string)
-importFrom(ggplot2,geom_text)
-importFrom(ggplot2,element_blank)
-importFrom(ggplot2,element_text)
-importFrom(ggplot2,scale_color_manual)
-importFrom(ggplot2,scale_linetype_manual)
-importFrom(ggplot2,scale_x_continuous)
-importFrom(ggplot2,theme_minimal)
-importFrom(ggplot2,scale_fill_viridis_c)
+import(ggplot2)
 importFrom(ggpubr,ggtexttable)
 importFrom(ggpubr,ggarrange)
 importFrom(ggpubr,annotate_figure)
@@ -182,6 +153,7 @@ importFrom(lhs,randomLHS)
 importFrom(magrittr,"%>%")
 importFrom(raster,crs)
 importFrom(raster,extract)
+importFrom(raster,extent)
 importFrom(raster,stack)
 importFrom(raster,mask)
 importFrom(raster,getValues)
@@ -223,4 +195,4 @@ importFrom(utils,combn)
 importFrom(utils,capture.output)
 importFrom(utils,head)
 importFrom(utils, install.packages)
-importFrom(utils,read.csv)
+importFrom(utils,read.csv)

NEWS.Rmd

@@ -0,0 +1,35 @@
+---
+title: "NEWS"
+author: "Dan Warren"
+date: "8/25/2020"
+output:
+  html_document:
+    df_print: paged
+---
+
+# ENMTools 1.x
+
+## ENMTools 1.0.1
+
+### Enhancements
+
+* Added variable importance tests via interface with the *vip* package
+* Added clamping for the *predict* functions, including plots of where clamping is happening
+* Added clamping for model construction functions, with a TRUE/FALSE switch defaulting to TRUE
+* Changed naming conventions for *predict* functions so that the suitability raster is in the $suitability slot, just as with modeling functions
+* Added progress bars for a lot of tests
+* Added "verbose" option for a lot of functions, defaulting to FALSE
+
+### Bug fixes
+
+* Added some code to fix a bug with recalibration.  Basically CalibratR's call to parallel was tanking on Mac OS, and there's just a little code snippet that needed to be added to the call to make that work.
+* Fixed interactive.plot generic and moved the function to its own file to make it easier to extend
+* Temporarily suppressing some warnings coming out of *leaflet* that are being produced by the recent *rgdal* changes
+* Fixed background sampling code to resample when necessary
+* Changed *enmtools.ranger* demo code to actually use *ranger* instead of *rf*
+* Fixed code for calculating p values for some of the hypothesis tests, the old code was getting wrong answers when there were repeated values
+
+
+## ENMTools 1.0.0
+
+* Initial CRAN release

NEWS.md

@@ -0,0 +1,40 @@
+ENMTools 1.x
+============
+
+ENMTools 1.0.1
+--------------
+
+### Enhancements
+
+-   Added variable importance tests via interface with the *vip* package
+-   Added clamping for the *predict* functions, including plots of where
+    clamping is happening
+-   Added clamping for model construction functions, with a TRUE/FALSE
+    switch defaulting to TRUE
+-   Changed naming conventions for *predict* functions so that the
+    suitability raster is in the $suitability slot, just as with
+    modeling functions
+-   Added progress bars for a lot of tests
+-   Added “verbose” option for a lot of functions, defaulting to FALSE
+
+### Bug fixes
+
+-   Added some code to fix a bug with recalibration. Basically
+    CalibratR’s call to parallel was tanking on Mac OS, and there’s just
+    a little code snippet that needed to be added to the call to make
+    that work.
+-   Fixed interactive.plot generic and moved the function to its own
+    file to make it easier to extend
+-   Temporarily suppressing some warnings coming out of *leaflet* that
+    are being produced by the recent *rgdal* changes
+-   Fixed background sampling code to resample when necessary
+-   Changed *enmtools.ranger* demo code to actually use *ranger* instead
+    of *rf*
+-   Fixed code for calculating p values for some of the hypothesis
+    tests, the old code was getting wrong answers when there were
+    repeated values
+
+ENMTools 1.0.0
+--------------
+
+-   Initial CRAN release

R/background.test.R

@@ -15,6 +15,8 @@
 #' @param bg.source Source for drawing background points.  If "points", it just uses the background points that are already in the species object.  If "range", it uses the range raster.  If "env", it draws points at randome from the entire study area outlined by the first environmental layer.
 #' @param low.memory When set to TRUE, replicate models are written to disc instead of being stored in the output object.  Replicate models stored in the output object contain paths to the replicate models on disk instead of the rasters themselves.
 #' @param rep.dir Directory for storing replicate models when low.memory is set to TRUE.  If not specified, the working directory will be used.
+#' @param verbose Controls printing of various messages progress reports.  Defaults to FALSE.
+#' @param clamp Controls whether empirical and replicate models should be clamped to the environment space used for training.
 #' @param ... Additional arguments to be passed to model fitting functions.
 #'
 #'
@@ -34,7 +36,7 @@
 #' f = pres ~ bio1 + bio12, nreps = 10)
 #' }
 
-background.test <- function(species.1, species.2, env, type, f = NULL, nreps = 99, test.type = "asymmetric", nback = 1000, bg.source = "default", low.memory = FALSE, rep.dir = NA, ...){
+background.test <- function(species.1, species.2, env, type, f = NULL, nreps = 99, test.type = "asymmetric", nback = 1000, bg.source = "default", low.memory = FALSE, rep.dir = NA, verbose = FALSE, clamp = TRUE, ...){
 
   # Build a description of the analysis to use for summaries and plot titles
   if(test.type == "symmetric"){
@@ -44,8 +46,8 @@ background.test <- function(species.1, species.2, env, type, f = NULL, nreps = 9
   }
   message(paste("\n", description, "\n"))
 
-  species.1 <- check.bg(species.1, env, nback = nback, bg.source = bg.source)
-  species.2 <- check.bg(species.2, env, nback = nback, bg.source = bg.source)
+  species.1 <- check.bg(species.1, env, nback = nback, bg.source = bg.source, verbose = verbose)
+  species.2 <- check.bg(species.2, env, nback = nback, bg.source = bg.source, verbose = verbose)
 
   # Check to make sure everything's okay
   background.precheck(species.1, species.2, env, type, f, nreps, test.type)
@@ -71,52 +73,75 @@ background.test <- function(species.1, species.2, env, type, f = NULL, nreps = 9
 
   # For starters we need to combine species background points so that each model
   # is being built with the same background
-  species.1$background.points <- rbind(species.1$background.points, species.2$background.points)
-  species.2$background.points <- rbind(species.1$background.points, species.2$background.points)
+  combined.background.points <- rbind(species.1$background.points, species.2$background.points)
+  species.1$background.points <- combined.background.points
+  species.2$background.points <- combined.background.points
+
+  # Clamping layers here so it's not done separately for every replicate
+  # and setting replicate clmaping to FALSE
+  if(clamp == TRUE){
+    combined.all.points <- rbind(species.1$presence.points, species.2$presence.points, combined.background.points)
+
+    # Adding env (skipped for BC otherwise)
+    this.df <- as.data.frame(extract(env, combined.all.points))
+
+    env <- clamp.env(this.df, env)
+  }
 
   # Building models for empirical data
   message("\nBuilding empirical models...\n")
   if(type == "glm"){
-    empirical.species.1.model <- enmtools.glm(species.1, env, f, ...)
-    empirical.species.2.model <- enmtools.glm(species.2, env, f, ...)
+    empirical.species.1.model <- enmtools.glm(species.1, env, f, clamp = FALSE, ...)
+    empirical.species.2.model <- enmtools.glm(species.2, env, f, clamp = FALSE, ...)
   }
 
   if(type == "gam"){
-    empirical.species.1.model <- enmtools.gam(species.1, env, f, ...)
-    empirical.species.2.model <- enmtools.gam(species.2, env, f, ...)
+    empirical.species.1.model <- enmtools.gam(species.1, env, f, clamp = FALSE, ...)
+    empirical.species.2.model <- enmtools.gam(species.2, env, f, clamp = FALSE, ...)
   }
 
   if(type == "mx"){
-    empirical.species.1.model <- enmtools.maxent(species.1, env, ...)
-    empirical.species.2.model <- enmtools.maxent(species.2, env, ...)
+    empirical.species.1.model <- enmtools.maxent(species.1, env, clamp = FALSE, ...)
+    empirical.species.2.model <- enmtools.maxent(species.2, env, clamp = FALSE, ...)
   }
 
   if(type == "bc"){
-    empirical.species.1.model <- enmtools.bc(species.1, env, ...)
-    empirical.species.2.model <- enmtools.bc(species.2, env, ...)
+    empirical.species.1.model <- enmtools.bc(species.1, env, clamp = FALSE, ...)
+    empirical.species.2.model <- enmtools.bc(species.2, env, clamp = FALSE, ...)
   }
 
   if(type == "dm"){
-    empirical.species.1.model <- enmtools.dm(species.1, env, ...)
-    empirical.species.2.model <- enmtools.dm(species.2, env, ...)
+    empirical.species.1.model <- enmtools.dm(species.1, env, clamp = FALSE, ...)
+    empirical.species.2.model <- enmtools.dm(species.2, env, clamp = FALSE, ...)
   }
 
   if(type == "rf"){
-    empirical.species.1.model <- enmtools.rf(species.1, env, ...)
-    empirical.species.2.model <- enmtools.rf(species.2, env, ...)
+    empirical.species.1.model <- enmtools.rf(species.1, env, clamp = FALSE, ...)
+    empirical.species.2.model <- enmtools.rf(species.2, env, clamp = FALSE, ...)
   }
 
   empirical.overlap <- c(unlist(raster.overlap(empirical.species.1.model, empirical.species.2.model)),
                          unlist(env.overlap(empirical.species.1.model, empirical.species.2.model, env = env)[1:3]))
   reps.overlap <- empirical.overlap
 
   message("\nBuilding replicate models...\n")
+
+  if (requireNamespace("progress", quietly = TRUE)) {
+    pb <- progress::progress_bar$new(
+      format = " [:bar] :percent eta: :eta",
+      total = nreps, clear = FALSE, width= 60)
+  }
+
   for(i in 1:nreps){
-    message(paste("\nReplicate", i, "...\n"))
+    if(verbose == TRUE){message(paste("\nReplicate", i, "...\n"))}
 
     rep.species.1 <- species.1
     rep.species.2 <- species.2
 
+    if (requireNamespace("progress", quietly = TRUE)) {
+      pb$tick()
+    }
+
     if(test.type == "symmetric"){
 
       #background sampling for species 1, only run when we're doing a symmetric test
@@ -137,33 +162,33 @@ background.test <- function(species.1, species.2, env, type, f = NULL, nreps = 9
 
     # Building the models for reps
     if(type == "glm"){
-      rep.species.1.model <- enmtools.glm(rep.species.1, env, f, ...)
-      rep.species.2.model <- enmtools.glm(rep.species.2, env, f, ...)
+      rep.species.1.model <- enmtools.glm(rep.species.1, env, f, clamp = FALSE, ...)
+      rep.species.2.model <- enmtools.glm(rep.species.2, env, f, clamp = FALSE, ...)
     }
 
     if(type == "gam"){
-      rep.species.1.model <- enmtools.gam(rep.species.1, env, f, ...)
-      rep.species.2.model <- enmtools.gam(rep.species.2, env, f, ...)
+      rep.species.1.model <- enmtools.gam(rep.species.1, env, f, clamp = FALSE, ...)
+      rep.species.2.model <- enmtools.gam(rep.species.2, env, f, clamp = FALSE, ...)
     }
 
     if(type == "mx"){
-      rep.species.1.model <- enmtools.maxent(rep.species.1, env, ...)
-      rep.species.2.model <- enmtools.maxent(rep.species.2, env, ...)
+      rep.species.1.model <- enmtools.maxent(rep.species.1, env, clamp = FALSE, ...)
+      rep.species.2.model <- enmtools.maxent(rep.species.2, env, clamp = FALSE, ...)
     }
 
     if(type == "bc"){
-      rep.species.1.model <- enmtools.bc(rep.species.1, env, ...)
-      rep.species.2.model <- enmtools.bc(rep.species.2, env, ...)
+      rep.species.1.model <- enmtools.bc(rep.species.1, env, clamp = FALSE, ...)
+      rep.species.2.model <- enmtools.bc(rep.species.2, env, clamp = FALSE, ...)
     }
 
     if(type == "dm"){
-      rep.species.1.model <- enmtools.dm(rep.species.1, env, ...)
-      rep.species.2.model <- enmtools.dm(rep.species.2, env, ...)
+      rep.species.1.model <- enmtools.dm(rep.species.1, env, clamp = FALSE, ...)
+      rep.species.2.model <- enmtools.dm(rep.species.2, env, clamp = FALSE, ...)
     }
 
     if(type == "rf"){
-      rep.species.1.model <- enmtools.rf(rep.species.1, env, ...)
-      rep.species.2.model <- enmtools.rf(rep.species.2, env, ...)
+      rep.species.1.model <- enmtools.rf(rep.species.1, env, clamp = FALSE, ...)
+      rep.species.2.model <- enmtools.rf(rep.species.2, env, clamp = FALSE, ...)
     }
 
     # Appending overlap to results
@@ -187,7 +212,7 @@ background.test <- function(species.1, species.2, env, type, f = NULL, nreps = 9
 
   rownames(reps.overlap) <- c("empirical", paste("rep", 1:nreps))
   print(reps.overlap)
-  p.values <- apply(reps.overlap, 2, function(x) 2 * (1 - max(mean(x > x[1]), mean(x < x[1]))))
+  p.values <- apply(reps.overlap, 2, function(x) min(rank(x)[1], rank(-x)[1])/length(x))
 
 
   d.plot <- qplot(reps.overlap[2:nrow(reps.overlap),"D"], geom = "histogram", fill = "density", alpha = 0.5) +