densitymodelling/dsmextra

Extrapolation assessments in density surface models

dsmextra's Issues

2 bugs in dsmextra?

Hi and Happy New Year!

Thanks for the great package - it's really helpful! I had started to cobble something much cruder together and so I was delighted when I found dsmextra!

I think I have found a couple of simple bugs. Apologies for not providing a complete reprex below, but I wasn't able to quickly pare down my segdata (327 rows) and predgrid (9124 rows) to make a minimal example that still displayed the behaviour in question.

1. Apply() sometimes returns the wrong type of structure

There are cases where apply() can return either a vector or matrix even though it was passed a dataframe, which messes up subsequent code.

There are two places in summarise_extrapolation.R where apply() is used to convert each column of a dataframe to character:

1. At line 144: resdf <- apply(resdf, 2, function(x) as.character(x)), and
2. At line 327: mic_resdf <- apply(mic_resdf, 2, function(x) as.character(x))

I have an example use case where apply() returns either a character vector or a 2-D character matrix on one of these rows instead of a dataframe, which messes up subsequent code that expects resdf or mic_resdf to be a dataframe.

For example, executing the following code where 'sst.sc' has no univariate or combinatorial extrapolation:

x <- compute_extrapolation(segments = segdata,
                           covariate.names = c("x.sc", "sst.sc"),
                           prediction.grid = predgrid,
                           coordinate.system = CRS("+proj=lcc +lon_0=-48 +lat_1=45.5 +lat_2=47 +ellps=GRS80"),
                           print.summary = TRUE,
                           save.summary = FALSE,
                           print.precision = 2)

causes the following output:

x <- compute_extrapolation(segments = segdata,
+                            covariate.names = c("x.sc", "sst.sc"),
+                            prediction.grid = predgrid,
+                            coordinate.system = CRS("+proj=lcc +lon_0=-48 +lat_1=45.5 +lat_2=47 +ellps=GRS80"),
+                            print.summary = TRUE,
+                            save.summary = FALSE,
+                            print.precision = 2)
Computing ...
Done!
 Show Traceback
 
 Rerun with Debug
 Error in rep("-----------", ncol(mic_resdf)) : invalid 'times' argument 

In this particular case, it fails at line 329 because apply() (at line 327) has returned a character vector and thus mic_resdf is no longer a dataframe at line 329.

Changing these two lines to use purrr::map_dfr() instead of apply() fixed the problem for me:

line 144: resdf <- purrr::map_dfr(resdf, as.character) and
line 327: mic_resdf <- purrr::map_dfr(mic_resdf, as.character)

2. max.univariate misspelled as max.univariates

There are four places in compare_covariates.R (on lines 245 and 250) where max.univariates is used, but I think you intended it to be max.univariate?

Line 245: min.univariate <- c(min.univariate, rep("-", times = length(max.univariates)-length(min.univariate)))

Line 250: max.univariates <- c(max.univariates, rep("-", times = length(min.univariate)-length(max.univariates)))

In my use case, executing the following:

compare_covariates(extrapolation.type = "both",
                   segments = segdata,
                   n.covariates = NULL,
                   covariate.names = c("x.sc", "sst.sc"),
                   prediction.grid = predgrid,
                   coordinate.system = CRS("+proj=lcc +lon_0=-48 +lat_1=45.5 +lat_2=47 +ellps=GRS80"),
                   create.plots = TRUE,
                   display.percent = TRUE)

causes the following error:

Computing ...
|=============================================================================================================|100% ~0 s remaining     

Creating summaries ...
Error in compare_covariates(extrapolation.type = "both", segments = segdata,  : 
  object 'max.univariates' not found

Changing max.univariates to max.univariate on these two lines fixed the problem for me.

Thanks!
Dave

Prediction grid cells are not regularly spaced

Some users have reported the below error when using compute_extrapolation (which also affects compute_nearby and extrapolation_analysis):

Error in compute_extrapolation(segments = Input, covariate.names = cov_vec2, : Prediction grid cells are not regularly spaced. A target raster resolution must be specified.

Typo in check_crs

There's a spelling error at line 263 of utils.R:

}else{supressWarnings(coordinate.system <- sp::CRS(coordinate.system))}

where "supressWarnings" is missing the second "p".

Thanks so much!
Dave

Package dependency issue

Hi Phil,

Len suggested looking at your cool new package for checking extrapolation. Unfortunately I have stumbled at the first hurdle. The packages WhatIf and Zelig were removed from CRAN on the 7th and 23rd March. I have found the archive versions and installed those but they need a whole sweep of dependencies.

For info Zelig needed the following packages and had to be installed prior to WhatIf:
install.packages(c("AER", "Amelia", "coda", "Formula", "geepack", "MatchIt", "maxLik", "MCMCpack", "survey", "VGAM"))

WhatIf required:
install.packages(c("lpSolve", "pbmcapply"))

Also, I'm using
R version 3.6.3 (2020-02-29)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows Server >= 2012 x64 (build 9200)

I'm up and running now though so can start to play.

Cheers,
Lindesay
(CREEM)

Function `safe_raster()` in `compute_extrapolation` not found

Hi there, I'm doing some trouble-shooting with my data (checking some unexpected results) so I've been using the source code for the compute_extrapolation and ExDet functions. However the function safe_raster from the compute_extraplation function is not available anywhere that I've looked. It looks to be a function just to generate rasters from the output list for univariate, analogue, and combinatorial layers but I wanted to highlight this just in case others want to use the source code.

compare_covariates, Error: Column names `extrapolation.univariate.n`...must not be duplicated

Running compare_covariates stops with the following error for both the spermwhales dataset or own data:

Error: Column names extrapolation.univariate.n, extrapolation.univariate.n, extrapolation.univariate.n, extrapolation.univariate.n, extrapolation.univariate.n, and 10 more must not be duplicated.
Run rlang::last_error() to see where the error occurred.

Output of rlang::last_error() attached below. I could not trace the issue based on the provided traceback.

Running the code elemnets the error is thrown in line 446 & 449

`if(extrapolation.type=="both"){

  var.ind.univariate <- purrr::map_df(.x = var.ind, .f = ~vars.univariate[.]) %>%
    tidyr::gather(data = .)
  
  var.ind.combinatorial <- purrr::map_df(.x = var.ind, .f = ~vars.combinatorial[.]) %>%
    tidyr::gather(data = .)`

I believe .f is expecting a numeric class without attr., and indeed the following works:
var.ind.univariate <- purrr::map_df(.x = var.ind, .f = ~as.numeric(vars.univariate[.])) %>% tidyr::gather(data = .)

Tracing the origin of the vars.univariate variable, it is based on the output of retrieve numbers of cells from the extrapolation, lines 213-219, with the following code:

vars <- exsum %>% purrr::map(.x = ., .f = ~.x[names(.x)%in%paste0(extype, ".n")]) %>% unlist(.)

Setting the function unlist(.), with unlist(. , use.names = FALSE) solves the issue for me.

sessionInfo()
R version 3.6.1 (2019-07-05)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 18362)

Matrix products: default

locale:
[1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United States.1252
[3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C
[5] LC_TIME=English_United States.1252

attached base packages:
[1] stats graphics grDevices utils datasets methods base

other attached packages:
[1] ggplot2_3.3.1 dsmextra_1.1.0 dismo_1.1-4 raster_3.1-5
[5] sp_1.4-2 sf_0.9-3 dplyr_1.0.0 ProjectTemplate_0.9.2

loaded via a namespace (and not attached):
[1] xfun_0.14 tidyselect_1.1.0 remotes_2.1.1 purrr_0.3.4 lattice_0.20-38
[6] colorspace_1.4-1 vctrs_0.3.1 generics_0.0.2 htmltools_0.4.0 utf8_1.1.4
[11] rlang_0.4.6 pkgbuild_1.0.8 startup_0.14.1 e1071_1.7-3 pillar_1.4.4
[16] glue_1.4.1 withr_2.2.0 DBI_1.1.0 lifecycle_0.2.0 plyr_1.8.6
[21] munsell_0.5.0 gtable_0.3.0 htmlwidgets_1.5.1 codetools_0.2-16 labeling_0.3
[26] knitr_1.28 callr_3.4.3 ps_1.3.3 crosstalk_1.1.0.1 curl_4.3
[31] class_7.3-15 fansi_0.4.1 highr_0.8 Rcpp_1.0.4.6 KernSmooth_2.23-15
[36] backports_1.1.7 scales_1.1.1 classInt_0.4-3 leaflet_2.0.3 farver_2.0.3
[41] packrat_0.4.8-1 digest_0.6.25 processx_3.4.2 cowplot_1.0.0 grid_3.6.1
[46] rprojroot_1.3-2 here_0.1 rgdal_1.5-8 cli_2.0.2 tools_3.6.1
[51] magrittr_1.5 tibble_3.0.1 crayon_1.3.4 tidyr_1.1.0 pkgconfig_2.0.3
[56] ellipsis_0.3.1 prettyunits_1.1.1 assertthat_0.2.1 rstudioapi_0.11 R6_2.4.1
[61] units_0.6-6 compiler_3.6.1

Running rlang::last_error() gives
22. stop(fallback)
21. signal_abort(cnd)
20. cnd_signal(cnd)
19. value[3L]
18. tryCatchOne(tryCatchList(expr, names[-nh], parentenv, handlers[-nh]), names[nh], parentenv, handlers[[nh]])
17. tryCatchList(expr, classes, parentenv, handlers)
16. tryCatch(force(expr), vctrs_error_names_cannot_be_empty = function(cnd) {cnd <- error_column_names_cannot_be_empty(detect_empty_names(name), parent = cnd)
cnd$body <- details ...
15. subclass_name_repair_errors(name = name, details = details, vec_as_names(name,
repair = .name_repair, quiet = quiet || !is_character(.name_repair)))
14.repaired_names(names2(x), .name_repair = .name_repair, quiet = quiet)
13.set_names(x, repaired_names(names2(x), .name_repair = .name_repair, quiet = quiet))
12.set_repaired_names(output, .name_repair = .name_repair)
11.tibble_quos(xs[!is_null], .rows, .name_repair)
10.tibble(!!!.x)
9.FUN(X[[i]], ...)
8.lapply(.x, .f, ...)
7.map(dots, function(.x) if (is.data.frame(.x)) .x else tibble(!!!.x))
6.dplyr::bind_rows(res, .id = .id)
5.purrr::map_df(.x = var.ind, .f = ~vars.univariate[.])
4.eval(lhs, parent, parent)
3.eval(lhs, parent, parent)
2.purrr::map_df(.x = var.ind, .f = ~vars.univariate[.]) %>% tidyr::gather(data = .)
1.compare_covariates(extrapolation.type = "both", covariate.names = my_cov,
n.covariates = NULL, segments = segs, prediction.grid = predgrid,
coordinate.system = my_crs, create.plots = TRUE, display.percent = TRUE)

compare_covariates, param requiring update

Function: compare_covariates

@param extrapolation.type Character string. Type of extrapolation to be assessed. Can be one of \code{univariate}, \code{combinatorial}, or \code{both} (default).

Shouldn't this read:
"multivariate" rather than combinatorial, as per lines 71-72

if(!extrapolation.type%in%c("both", "univariate", "multivariate"))
stop("Unknown extrapolation type")

proj_rasters fails if univariate (or analogue) are NULL

Hi Phil,

I discovered that proj_rasters() fails if either ll$univariate or ll$analogue are NULL. I have a dataset where there is no univariate extrapolation so ll$univariate will be NULL. On another occasion, I messed up the dataframes passed to extrapolation_analysis() and ended up with no analogue either, and proj_rasters() failed in this case too.

As a workaround in my local copy of the package, I have passed the "types" object from map_extrapolation() to proj_rasters ().

Then within proj_rasters() I wrap two parts of the code in a conditional to test the value of types. Not sure if this is sensible or not, but it seems to fix the problem.

proj_rasters <- function(ll, coordinate.system, types){

  suppressWarnings(crs.webmerc <- sp::CRS("+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext  +no_defs"))

  llr <- ll # Copy list

  if(all(c("univariate", "analogue") %in% types)) {
      # Univariate extrapolation is negative by definition
      # When only a small number of cells are subject to UE, the resampling
      # may result in the loss of some of them.
      # By recording the indices of UE cells, we can perform a simplistic
      # correction to make sure they show up on the map.
    
      analogue.xy <- raster::as.data.frame(llr$analogue, xy = TRUE) %>% stats::na.omit(.)
      analogue.xy <- sp::SpatialPointsDataFrame(coords = analogue.xy[, c("x", "y")],
                                                data = analogue.xy,
                                                proj4string = coordinate.system)
      analogue.xy <- sp::spTransform(analogue.xy, CRSobj = crs.webmerc)
    
      univariate.ind <- raster::Which(llr$univariate < 0, cells = TRUE)
      univariate.values <- llr$univariate[univariate.ind]
    
      univariate.xy <- raster::as.data.frame(llr$univariate, xy = TRUE) %>% stats::na.omit(.)
      univariate.xy <- sp::SpatialPointsDataFrame(coords = univariate.xy[, c("x", "y")],
                                                  data = univariate.xy,
                                                  proj4string = coordinate.system)
      univariate.xy <- sp::spTransform(univariate.xy, CRSobj = crs.webmerc)
    }
    
    llr$all <- NULL
    llr <- purrr::discard(.x = llr, is.null)
  
    suppressWarnings(
    llr <- purrr::map(.x = llr, # Same extent as the full raster, allows correct alignment
                      .f = ~raster::projectRaster(from = .x,
                                                  to = ll$all,
                                                  method = 'ngb')) %>%
      purrr::map(.x = ., # CRS used by leaflet
                 .f = ~raster::projectRaster(from = .,
                                             crs = crs.webmerc,
                                             method = 'ngb')))

  if(all(c("univariate", "analogue") %in% types)) {
    llr.univariate.ind <- raster::cellFromXY(object = llr$univariate,
                                             xy = sp::coordinates(univariate.xy))
  
    llr$univariate[llr.univariate.ind[which(is.na(llr$univariate[llr.univariate.ind]))]] <-    univariate.values[which(is.na(llr$univariate[llr.univariate.ind]))]
  
    r1 <- raster::as.data.frame(llr$univariate, xy = TRUE)
    r2 <- raster::as.data.frame(llr$analogue, xy = TRUE)
    names(r1) <- names(r2) <- c("x", "y", "ExDet")
  
    duplicate.cells <- rbind(r1, r2) %>%
      stats::na.omit(.) %>%
      dplyr::select(., x, y) %>%
      .[duplicated(.),]
  
    llr.analogue.ind <- raster::cellFromXY(object = llr$analogue,
                                           xy = duplicate.cells)
  
  
    llr$analogue[llr.analogue.ind] <- NA
  }

  return(llr)}

Best!
Dave

ExDet function not found

Hi,
I am very interested in the function "ExDet" that you developed. I used these lines to install "dsmextra" package :

if (!require("remotes")) install.packages("remotes")
remotes::install_github("densitymodelling/dsmextra", force = TRUE)

I installed the latest version of each additional package needed and I ran library(dsmextra).

The package is well-installed but when I run "ExDet" function I got : ExDet(ref, tg, xp) : could not find function "ExDet"
Any idea to fix my issue ?

Many thanks in advance
Jade