NA value for stratum area in simulated data

[LHMarshall]

reproducible example (with appropriate shapefile)

Seems to be related to an object being detected in one stratum from a sampler in a different stratum.

# load simulation package
library(dsims)
library(sf)

#Load shapefile Of the region
shapefile.name <-st_read(dsn="DUStrataSeptember.shp")

#Make sure the projection is in equal area if not dssd become unstable.
sf::st_crs(shapefile.name)
# Define WGS 84 / North Pole LAEA Canada (EPSG: 3573)
proj4string<-"+proj=laea +lat_0=90 +lon_0=-100 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs "
# Project the study area on to a flat plane
projected.strata <- sf::st_transform(shapefile.name, crs = proj4string)

#Make Region
region <- make.region(region.name = "Study.Area",
                      shape = projected.strata)
# There are 15 strata you had only provided 14 names
region <- make.region(region.name = "Study.Area",
                      strata.name = c("A","B","C","D","E","F","G","H","I","J","K","L","M","N", "O"),
                      shape = projected.strata)
##Maybe the mistake is the strata.name? Dsims just overwrite my "Name" column in my shapefile. 
plot(region)

#Creating a Coverage Grid
sf_use_s2(FALSE)
library(lwgeom)
cover<-make.coverage(region, n.grid.points = 1000)
plot(region, cover)

#Design transect
default.design <- make.design(region = region,
                              transect.type = "line",
                              design = "systematic",
                              spacing = c(13500,5500,10000,14000,8500,8500,4000,9500,14000,14000,8000,12000,15000,15000,8000),
                              design.angle = c(05,35,75,75,45,165,0,45,75,135,45,155,10,155,45),
                              edge.protocol = "minus",
                              truncation = 1000, 
                              coverage.grid = cover)
transects <- generate.transects(default.design)
plot(region, transects)

density <- make.density(region = region, 
                        x.space = 5000, 
                        constant = c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1))

#Area of high/low density are based on point density analysis ran on groups 
# Add some areas of higher 
density <- add.hotspot(object = density,
                       centre = c(-442000, -2330000),
                       sigma = 30000,
                       amplitude = 10)
# Add some areas of higher 
density <- add.hotspot(object = density,
                       centre = c(-360000, -2320000),
                       sigma = 10000,
                       amplitude = 10)
# Add some areas of higher 
density <- add.hotspot(object = density,
                       centre = c(-300000, -2300000),
                       sigma = 10000,
                       amplitude = 10)

# Plot the density surface
plot(density,region)


# Create the population description, with a population size N = 3,825 and need to be listed by strata; fixed.N = True the population
#is generated from the values of N otherwise it is generated from the density grid.
pop.desc <- make.population.description(region = region, 
                                        density = density,
                                        N = c(1,821,111,400,470,89,487,5,0,65,130,236,0,0,10),
                                        fixed.N = TRUE)

#Detectability
# Make a simple half normal detection function with a scale parameter of 225
detect.hn <- make.detectability(key.function = "hn",
                                scale.param = 350, 
                                truncation = 1000)
# We can now visualise these detection functions
plot(detect.hn, pop.desc)

#Analysis
ddf.analyses <- make.ds.analysis(dfmodel = list(~1),
                                 key = c("hn", "hr"),
                                 criteria = "AIC",
                                 truncation = 1000)

#Simulation
sim <- make.simulation(reps = 2, 
                       design = default.design,
                       population.description = pop.desc,
                       detectability = detect.hn,
                       ds.analysis = ddf.analyses)


set.seed(136)
survey <- run.survey(sim)
index <- which(is.na(survey@dist.data$Area))
index
survey@dist.data[325:330,]

#    object individual Region.Label Sample.Label  distance         x        y
#274    274       2230            G          144  19.80172 -414361.7 -2329621
#275    275        164            B          145  33.44826 -410414.9 -2303713
#276    276       1985            G          145  69.96048 -410451.5 -2327677
#277    277       2258            G          145 116.64116 -410264.9 -2318848
#278    278       2317            G          145 165.36364 -410216.1 -2310655
#279    279       2354            G          145 505.76041 -410887.3 -2334331
#          Effort       Area
#274 32059.96 [m] 7032731502
#275       NA [m]         NA
#276 33209.73 [m] 7032731502
#277 33209.73 [m] 7032731502
#278 33209.73 [m] 7032731502
#279 33209.73 [m] 7032731502

[LHMarshall]

Now ~70% of the time there are no detections? Likely either related to the fix of this bug or #80