diff --git a/DESCRIPTION b/DESCRIPTION index d04e625..f6b9143 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -21,10 +21,9 @@ Suggests: knitr, lwgeom, rmarkdown -VignetteBuilder: knitr Type: Package Title: Distance Sampling Simulations -Version: 1.0.4.9000 +Version: 1.0.5 Authors@R: c( person("Laura", "Marshall", email = "lhm@st-and.ac.uk", role = c("aut", "cre")), person("Thomas", "Len", email = "len.thomas@st-andrews.ac.uk", role = "ctb")) @@ -43,7 +42,7 @@ Description: Performs distance sampling simulations. 'dsims' repeatedly generate 978-0199225873). General distance sampling methods are detailed in Introduction to Distance Sampling: Estimating Abundance of Biological Populations, Buckland et. al. (2004, ISBN-13: 978-0198509271). Find out more about estimating - animal/plant abundance with distance sampling at . + animal/plant abundance with distance sampling at . License: GPL (>=2) Language: en-GB URL: https://github.com/DistanceDevelopment/dsims @@ -93,6 +92,7 @@ Collate: 'rztpois.R' 'save.sim.results.R' 'simulate.detections.R' + 'simulation.consistency.check.R' 'single.sim.loop.R' 'store.ddf.results.R' 'store.dht.results.R' diff --git a/NEWS.md b/NEWS.md index ca3410b..39062bf 100644 --- a/NEWS.md +++ b/NEWS.md @@ -3,6 +3,7 @@ Bug Fixes * Fixed bug which generated NA's as scale parameters when factor covariates were included. Issue #89 +* Simulation validation checks: consistency in truncation distances (Issue #76), consistency in region (Issue #88) # dsims 1.0.4 diff --git a/R/ClassConstructors.R b/R/ClassConstructors.R index faf055b..c669fbc 100755 --- a/R/ClassConstructors.R +++ b/R/ClassConstructors.R @@ -31,7 +31,7 @@ #' @importFrom dssd make.region #' @importFrom methods new #' @author Laura Marshall -#' @seealso \code{\link{make.region}} +#' @seealso \code{\link[dssd]{make.region}} #' @examples #' # A simple density surface with a constant value of 1 can be created within a rectangular #' # Create a region from shapefile @@ -112,7 +112,7 @@ make.density <- function(region = make.region(), x.space = 20, y.space = NULL, c #' lognormal \tab meanlog \tab sdlog \cr #' } #' -#' @param region the Region object in which this population exists (see \link{make.region}). +#' @param region the Region object in which this population exists (see \link[dssd]{make.region}). #' @param density the Density object describing the distribution of the individuals / clusters (see \link{make.density}). #' @param covariates Named list with one named entry per individual-level covariate. Cluster sizes can be defined here, it must be named 'size'. The distribution of covariate values can either be defined by specifying a particular distribution and its parameters or as a discrete distribution in a dataframe. Dataframes should have columns level and prob (and optionally strata) specifying the covariates levels, probabilities and strata if they are strata specific. Distributions can be defined as lists with named entries distribution and the relevant parameters as specified in details. A list of distributions can be provided with one for each strata. #' @param N the number of individuals / clusters in a population with one value per @@ -377,7 +377,7 @@ make.detectability <- function(key.function = "hn", scale.param = 25, shape.para #' @export #' @importFrom methods new is #' @author Laura Marshall -#' @seealso \code{\link{ds}} \code{\link{make.simulation}} +#' @seealso \code{\link[Distance]{ds}} \code{\link{make.simulation}} #' @examples #' #' # Model selection considering both a half-normal and a hazard-rate model @@ -465,7 +465,7 @@ make.ds.analysis <- function(dfmodel = list(~1), #' line) it can run a simple simulation example. See examples. #' @param reps number of times the simulation should be repeated #' @param design an object of class Survey.Design created by a call to -#' \link{make.design} +#' \link[dssd]{make.design} #' @param population.description an object of class Population.Description #' created by a call to \link{make.population.description} #' @param detectability and object of class Detectability created by a call to @@ -477,7 +477,7 @@ make.ds.analysis <- function(dfmodel = list(~1), #' @importFrom methods new is #' @importFrom dssd make.region make.design #' @author Laura Marshall -#' @seealso \code{\link{make.region}} \code{\link{make.density}} \code{\link{make.population.description}} \code{\link{make.detectability}} \code{\link{make.ds.analysis}} \code{\link{make.design}} +#' @seealso \code{\link[dssd]{make.region}} \code{\link{make.density}} \code{\link{make.population.description}} \code{\link{make.detectability}} \code{\link{make.ds.analysis}} \code{\link[dssd]{make.design}} #' @examples #' # Create a basic rectangular study area #' region <- make.region() diff --git a/R/Simulation.R b/R/Simulation.R index c65dd59..690df96 100755 --- a/R/Simulation.R +++ b/R/Simulation.R @@ -138,6 +138,8 @@ setMethod( .Object@ds.analysis <- ds.analysis .Object@results <- results .Object@warnings <- list() + #Consistency check + .Object <- simulation.consistency.check(.Object) #Check object is valid valid <- validObject(.Object, test = TRUE) if(is(valid, "character")){ @@ -150,15 +152,34 @@ setMethod( setValidity("Simulation", function(object){ + strata.names <- object@design@region@strata.name - # truncation - # Check to see if the analysis truncation distance is larger than the - # if(length(object@ds.analysis@truncation[[1]]) > 0){ - # if(object@ds.analysis@truncation > object@detectability@truncation || - # object@ds.analysis@truncation > object@design@truncation){ - # warning("The truncation distance for analysis is larger than the truncation distance for data generation, this will likely cause biased results.", immediate. = TRUE, call. = FALSE) - # } - # } + + # REGION CHECKS + # Check the intersection between design and population regions - they must overlap + intersect <- suppressWarnings(sf::st_intersection(object@design@region@region, object@population.description@density@density.surface[[1]])) + intersect.area <- sum(sf::st_area(intersect)) + if(intersect.area == 0){ + return("The regions associated with the design and the population description do not overlap!") + } + # Compare the area of the design region and population region with the intersection. + pop.area <- sum(sf::st_area(object@population.description@density@density.surface[[1]])) + design.area <- sum(sf::st_area(object@design@region@region)) + area.diff.pop <- abs(pop.area-intersect.area) + area.diff.des <- abs(design.area-intersect.area) + # Check that the areas of both the population and design and their intersection match, using a toleance of 0.1%. + if(area.diff.pop > 0.001*pop.area || area.diff.des > 0.001*design.area){ + return("The regions for the population density surface and the design do not match, these must be the same.") + } + # The following checks are for when simulations allow differing regions between population and design - this cannot be implemented yet as requires updates to population generation code + # Using 0.5% as toerance for difference + if(area.diff.pop > 0.005*pop.area && area.diff.des < 0.005*design.area){ + warning("The population density surface extends beyond the design survey region, only part of the population will be surveyed.", immediate. = TRUE, call. = FALSE) + } + if(intersect.area < design.area && area.diff.des > 0.005*design.area){ + warning("The survey design extends beyond the density grid for the population, some survey areas will have no animals.", immediate. = TRUE, call. = FALSE) + } + # Population.Description checks pop.desc <- object@population.description diff --git a/R/Survey.LT.R b/R/Survey.LT.R index af6f232..c5b50d4 100755 --- a/R/Survey.LT.R +++ b/R/Survey.LT.R @@ -12,7 +12,7 @@ #' maximum distance from the transect at which animals may be detected. #' @keywords classes #' @importClassesFrom dssd Line.Transect -#' @seealso \code{\link{make.design}} +#' @seealso \code{\link[dssd]{make.design}} #' @export setClass(Class = "Survey.LT", representation = representation(transect = "Line.Transect", diff --git a/R/Survey.PT.R b/R/Survey.PT.R index 9ba2ad5..2065753 100755 --- a/R/Survey.PT.R +++ b/R/Survey.PT.R @@ -12,7 +12,7 @@ #' maximum distance from the transect at which animals may be detected. #' @keywords classes #' @importClassesFrom dssd Point.Transect -#' @seealso \code{\link{make.design}} +#' @seealso \code{\link[dssd]{make.design}} #' @export setClass(Class = "Survey.PT", representation = representation(transect = "Point.Transect", diff --git a/R/dsims-package.R b/R/dsims-package.R index 534853e..5278318 100644 --- a/R/dsims-package.R +++ b/R/dsims-package.R @@ -18,11 +18,11 @@ #' package relied on survey transects already being contained in shapefiles within #' the supplied directory, dsims will generate the survey transects directly in R. #' -#' The main functions in this package are: \link{make.density}, \link{make.population.description}, \link{make.detectability}, \link{make.ds.analysis}, \link{make.simulation}, \link{run.survey} and \link{run.simulation}. See also \link{make.region} and \link{make.design} in the dssd package for examples of how to define study regions and designs. +#' The main functions in this package are: \link{make.density}, \link{make.population.description}, \link{make.detectability}, \link{make.ds.analysis}, \link{make.simulation}, \link{run.survey} and \link{run.simulation}. See also \link[dssd]{make.region} and \link[dssd]{make.design} in the dssd package for examples of how to define study regions and designs. #' -#' Further information on distance sampling methods and example code is available at \url{http://distancesampling.org/R/}. -#' -#' Also see our website for vignettes / example code at \url{http://examples.distancesampling.org}. +#' Further information on distance sampling methods and example code is +#' available at \url{https://distancesampling.org/}. Specifically, see our +#' website for vignettes / example code at \url{https://distancesampling.org/resources/vignettes.html}. #' #' For help with distance sampling and this package, there is a Google Group \url{https://groups.google.com/forum/#!forum/distance-sampling}. #' diff --git a/R/simulation.consistency.check.R b/R/simulation.consistency.check.R new file mode 100644 index 0000000..025ebcc --- /dev/null +++ b/R/simulation.consistency.check.R @@ -0,0 +1,24 @@ +simulation.consistency.check <- function(sim){ + # This function checks for consistency across objects within the simulation + + # TRUNCATION CHECKS + # The truncation distances should be the same for design and detectability + dd.trunc.match <- ifelse(sim@design@truncation == sim@detectability@truncation, TRUE, FALSE) + analysis.trunc.greater <- ifelse(sim@ds.analysis@truncation[[1]] > sim@detectability@truncation, TRUE, FALSE) + + # If there is a design/detect truncation mismatch but analysis truncation distance ok + if(!dd.trunc.match && !analysis.trunc.greater){ + warning(paste("Truncation distance for design and detectability differ, updating design truncation to be ", sim@detectability@truncation, ".", sep = ""), immediate. = TRUE, call. = FALSE) + sim@design@truncation <- sim@detectability@truncation + # If design/detect truncations match but analysis truncation is larger + }else if(dd.trunc.match && analysis.trunc.greater){ + warning("Truncation distance for analysis is larger than for design/detectability this may introduce bias!", immediate. = TRUE, call. = FALSE) + # If there is both design/detect mismatch and analysis truncation is larger + }else if(!dd.trunc.match && analysis.trunc.greater){ + warning(paste("Truncation distance for design and detectability differ, updating design truncation to be ", sim@detectability@truncation, ". In addition, analysis truncation is greater than ", sim@detectability@truncation, " this may introduce bias!", sep = ""), immediate. = TRUE, call. = FALSE) + sim@design@truncation <- sim@detectability@truncation + } + + # Return object with any updates + return(sim) +} \ No newline at end of file diff --git a/README.md b/README.md index 445dcd6..5ed7df1 100644 --- a/README.md +++ b/README.md @@ -20,7 +20,8 @@ There is currently three vignette within the dsims package to help you get start The easiest way to get `dsims` is to install it from CRAN within R-studio or the R interface. We endeavour to make all new functionality available on CRAN in a timely manor. However, if you wish to download the development version with the latest updates immediately you can do this using Hadley Wickham's `devtools` package: - install.packages("devtools") + # First, ensure you have a copy of the `devtools` package: + if (!nzchar(system.file(package = "devtools"))) install.packages("devtools") then install `dsims` from github: diff --git a/man/Survey.LT-class.Rd b/man/Survey.LT-class.Rd index cde60ca..787d4e2 100644 --- a/man/Survey.LT-class.Rd +++ b/man/Survey.LT-class.Rd @@ -19,6 +19,6 @@ maximum distance from the transect at which animals may be detected.} }} \seealso{ -\code{\link{make.design}} +\code{\link[dssd]{make.design}} } \keyword{classes} diff --git a/man/Survey.PT-class.Rd b/man/Survey.PT-class.Rd index 7648e34..e8fb63f 100644 --- a/man/Survey.PT-class.Rd +++ b/man/Survey.PT-class.Rd @@ -19,6 +19,6 @@ maximum distance from the transect at which animals may be detected.} }} \seealso{ -\code{\link{make.design}} +\code{\link[dssd]{make.design}} } \keyword{classes} diff --git a/man/dsims-package.Rd b/man/dsims-package.Rd index 6ce2a28..67938b6 100644 --- a/man/dsims-package.Rd +++ b/man/dsims-package.Rd @@ -24,11 +24,11 @@ regions, designs and generate the survey transects. While the 'DSsim' simulation package relied on survey transects already being contained in shapefiles within the supplied directory, dsims will generate the survey transects directly in R. - The main functions in this package are: \link{make.density}, \link{make.population.description}, \link{make.detectability}, \link{make.ds.analysis}, \link{make.simulation}, \link{run.survey} and \link{run.simulation}. See also \link{make.region} and \link{make.design} in the dssd package for examples of how to define study regions and designs. + The main functions in this package are: \link{make.density}, \link{make.population.description}, \link{make.detectability}, \link{make.ds.analysis}, \link{make.simulation}, \link{run.survey} and \link{run.simulation}. See also \link[dssd]{make.region} and \link[dssd]{make.design} in the dssd package for examples of how to define study regions and designs. -Further information on distance sampling methods and example code is available at \url{http://distancesampling.org/R/}. - -Also see our website for vignettes / example code at \url{http://examples.distancesampling.org}. +Further information on distance sampling methods and example code is +available at \url{https://distancesampling.org/}. Specifically, see our +website for vignettes / example code at \url{https://distancesampling.org/resources/vignettes.html}. For help with distance sampling and this package, there is a Google Group \url{https://groups.google.com/forum/#!forum/distance-sampling}. } diff --git a/man/make.density.Rd b/man/make.density.Rd index 17699d1..32cd623 100644 --- a/man/make.density.Rd +++ b/man/make.density.Rd @@ -76,7 +76,7 @@ plot(density, region) } \seealso{ -\code{\link{make.region}} +\code{\link[dssd]{make.region}} } \author{ Laura Marshall diff --git a/man/make.ds.analysis.Rd b/man/make.ds.analysis.Rd index 15ddfff..307adbf 100644 --- a/man/make.ds.analysis.Rd +++ b/man/make.ds.analysis.Rd @@ -95,7 +95,7 @@ ds.analyses <- make.ds.analysis(dfmodel = list(~1, ~size, ~1), } \seealso{ -\code{\link{ds}} \code{\link{make.simulation}} +\code{\link[Distance]{ds}} \code{\link{make.simulation}} } \author{ Laura Marshall diff --git a/man/make.population.description.Rd b/man/make.population.description.Rd index 2254064..ddbeaa4 100644 --- a/man/make.population.description.Rd +++ b/man/make.population.description.Rd @@ -13,7 +13,7 @@ make.population.description( ) } \arguments{ -\item{region}{the Region object in which this population exists (see \link{make.region}).} +\item{region}{the Region object in which this population exists (see \link[dssd]{make.region}).} \item{density}{the Density object describing the distribution of the individuals / clusters (see \link{make.density}).} diff --git a/man/make.simulation.Rd b/man/make.simulation.Rd index 77a0b29..fb6837d 100644 --- a/man/make.simulation.Rd +++ b/man/make.simulation.Rd @@ -16,7 +16,7 @@ make.simulation( \item{reps}{number of times the simulation should be repeated} \item{design}{an object of class Survey.Design created by a call to -\link{make.design}} +\link[dssd]{make.design}} \item{population.description}{an object of class Population.Description created by a call to \link{make.population.description}} @@ -105,7 +105,7 @@ vignette("GettingStarted", 'dsims') } \seealso{ -\code{\link{make.region}} \code{\link{make.density}} \code{\link{make.population.description}} \code{\link{make.detectability}} \code{\link{make.ds.analysis}} \code{\link{make.design}} +\code{\link[dssd]{make.region}} \code{\link{make.density}} \code{\link{make.population.description}} \code{\link{make.detectability}} \code{\link{make.ds.analysis}} \code{\link[dssd]{make.design}} } \author{ Laura Marshall diff --git a/tests/testthat/test-ProblemCases.R b/tests/testthat/test-ProblemCases.R index 790f281..d1047db 100644 --- a/tests/testthat/test-ProblemCases.R +++ b/tests/testthat/test-ProblemCases.R @@ -98,4 +98,45 @@ test_that("AICc simulation", { sim <- run.simulation(sim) expect_s4_class(sim, "Simulation") +}) + +test_that("Check cannot get detections greater than truncation", { + + outer <- matrix(c(0,0,15,0,15,10,0,10,0,0),ncol=2, byrow=TRUE) + pol1 <- sf::st_polygon(list(outer)) + pol2 <- sf::st_polygon(list(outer + 15)) + pol3 <- sf::st_polygon(list(outer + 30)) + sfc <- sf::st_sfc(pol1,pol2,pol3) + strata.names <- c("SW", "central", "NE") + mp1 <- sf::st_sf(strata = strata.names, geom = sfc) + + region <- make.region(region.name = "study.area", + strata.name = strata.names, + shape = mp1) + + density <- make.density(region = region, + x.space = 0.25, + constant = rep(1,3)) + + popdesc <- make.population.description(region = region, + density = density, + N = rep(100,3), + fixed.N = TRUE) + + design <- make.design(region, + spacing = 1, + truncation = 0.5) + + detect <- make.detectability(scale.param = 0.5, + truncation = 0.5) + + suppressWarnings(sim <- make.simulation(reps = 10, + design = design, + population.description = popdesc, + detectability = detect)) + + + survey <- run.survey(sim) + # Check all distances are less than truncation distance (Issue #76) + expect_true(all(survey@dist.data$distance <= 0.5)) }) \ No newline at end of file diff --git a/tests/testthat/test-check_Constructors.R b/tests/testthat/test-check_Constructors.R index 8815d04..bdfb2d1 100644 --- a/tests/testthat/test-check_Constructors.R +++ b/tests/testthat/test-check_Constructors.R @@ -276,3 +276,127 @@ test_that("Can create objects or return correct error / warning messages", { survey <- run.survey(sim) }) + +test_that("Test simulation validation and consistency checks", { + # Create a basic rectangular study area + region <- make.region() + + # Make a density grid (large spacing for speed) + density <- make.density(region = region, + x.space = 300, + y.space = 100, + constant = 1) + density <- add.hotspot(density, centre = c(1000, 100), sigma = 250, amplitude = 10) + + # Define the population description + popdsc <- make.population.description(region = region, + density = density, + N = 200) + + # Define the detecability + detect <- make.detectability(key.function = "hn", + scale.param = 25, + truncation = 50) + + # Define the design + design <- make.design(region = region, + transect.type = "line", + design = "systematic", + samplers = 20, + design.angle = 0, + truncation = 1) + + # Define the analyses + ds.analyses <- make.ds.analysis(dfmodel = ~1, + key = "hn", + truncation = 100, + criteria = "AIC") + + # Test when truncation design/detect mismatch and analysis truncation too large + expect_warning(simulation <- make.simulation(reps = 1, + design = design, + population.description = popdsc, + detectability = detect, + ds.analysis = ds.analyses), + "Truncation distance for design and detectability differ, updating design truncation to be 50. In addition, analysis truncation is greater than 50 this may introduce bias!") + + design@truncation <- 50 + # Test when only analysis truncation too large + expect_warning(simulation <- make.simulation(reps = 1, + design = design, + population.description = popdsc, + detectability = detect, + ds.analysis = ds.analyses), + "Truncation distance for analysis is larger than for design/detectability this may introduce bias!") + + design@truncation <- 20 + ds.analyses@truncation[[1]] <- 50 + # Test when truncation design/detect mismatch + expect_warning(simulation <- make.simulation(reps = 1, + design = design, + population.description = popdsc, + detectability = detect, + ds.analysis = ds.analyses), + "Truncation distance for design and detectability differ, updating design truncation to be 50.") + + # Test case where there is no overlap between region used for population and region used for design. + outer <- matrix(c(0,1000,2000,1000,2000,1500,0,1500,0,1000),ncol=2, byrow=TRUE) + mp1 <- sf::st_polygon(list(outer)) + region <- make.region(shape=mp1, units = "km") + pop.desc <- make.population.description(region = region, + density = make.density(region=region), + N=200) + + expect_error(make.simulation(population.description = pop.desc), + "The regions associated with the design and the population description do not overlap!") + + # Test case where the design only covers part of the population. + outer <- matrix(c(0,0,2000,0,2000,1500,0,1500,0,0),ncol=2, byrow=TRUE) + mp1 <- sf::st_polygon(list(outer)) + region <- make.region(shape=mp1, units = "km") + pop.desc <- make.population.description(region = region, + density = make.density(region=region), + N=200) + + expect_error(make.simulation(population.description = pop.desc), + "The regions for the population density surface and the design do not match, these must be the same.") + + #expect_warning(make.simulation(population.description = pop.desc),"The population density surface extends beyond the design survey region, only part of the population will be surveyed.") + + + # Test case where the design only covers part of the population. + outer <- matrix(c(-100,-100,2500,-100,2500,800,-100,800,-100,-100),ncol=2, byrow=TRUE) + mp1 <- sf::st_polygon(list(outer)) + region <- make.region(shape=mp1, units = "km") + pop.desc <- make.population.description(region = region, + density = make.density(region=region), + N=200) + expect_error(make.simulation(population.description = pop.desc), + "The regions for the population density surface and the design do not match, these must be the same.") + + #expect_warning(make.simulation(population.description = pop.desc),"The population density surface extends beyond the design survey region, only part of the population will be surveyed.") + + # Test case where the design extends beyond the population. + outer <- matrix(c(500,100,1500,100,1500,400,500,400,500,100),ncol=2, byrow=TRUE) + mp1 <- sf::st_polygon(list(outer)) + region <- make.region(shape=mp1, units = "km") + pop.desc <- make.population.description(region = region, + density = make.density(region=region), + N=200) + expect_error(make.simulation(population.description = pop.desc), + "The regions for the population density surface and the design do not match, these must be the same.") + #expect_warning(make.simulation(population.description = pop.desc),"The survey design extends beyond the density grid for the population, some survey areas will have no animals.") + + # Test case where the design extends beyond the population. + outer <- matrix(c(1500,100,2500,100,2500,400,1500,400,1500,100),ncol=2, byrow=TRUE) + mp1 <- sf::st_polygon(list(outer)) + region <- make.region(shape=mp1, units = "km") + pop.desc <- make.population.description(region = region, + density = make.density(region=region), + N=200) + + expect_error(make.simulation(population.description = pop.desc), + "The regions for the population density surface and the design do not match, these must be the same.") + + #expect_warning(make.simulation(population.description = pop.desc),"The population density surface extends beyond the design survey region, only part of the population will be surveyed.") +}) diff --git a/tests/testthat/test-check_ExtremeCases.R b/tests/testthat/test-check_ExtremeCases.R index 1a1673d..1287c78 100644 --- a/tests/testthat/test-check_ExtremeCases.R +++ b/tests/testthat/test-check_ExtremeCases.R @@ -14,7 +14,7 @@ test_that("Test very few or no detections", { N = 1) detect <- make.detectability(key.function = "hn", scale.param = 0.000000001, - truncation = 0.00001) + truncation = 0.000001) design <- make.design(region = region, transect.type = "line", samplers = 1,