In MurrayEfford/secrBVN: Bivariate Normal Home Range SECR
\vspace{16pt}
The small package secrBVN is used to evaluate the performance of SECR estimators when home ranges are bivariate normal (BVN) or uniform (flat-topped) ellipses. We assume detection hazard is directly proportional to home range utilisation (activity). Code to use secrBVN for the simulations of Efford (in prep.) is provided in the Appendix.
The key user-visible functions are simpopn.bvn, plotpopn.bvn, simcapt.bvn, runEllipseSim, simsum, simplot and anisotropic.fit.
\vspace{12pt}
Generating and plotting elliptical home ranges
simpopn.bvn is a wrapper for the secr function sim.popn that adds attributes specifying a bivariate normal home range shape, size and orientation for each individual. By default, shape and size are the same for all individuals, but a mechanism is provided to vary them individually (see Heterogeneous elliptical home ranges).
First, load the package.
library(secrBVN)
simfolder <- "d:/density communication/noncircularity/paper/simulations/"
runall <- FALSE # skip lengthy simulations
tempgrid <- make.grid(nx = 10, ny = 10)
par(mfrow = c(2,4), mar = c(2,2,2.6,2), xpd = TRUE)
for (i in 1:4) {
s2xy <- 25^2 * c(1/i, i)
random.pop <- simpopn.bvn(s2xy = s2xy, core = tempgrid, buffer = 100, D = 1)
plotpopn.bvn(random.pop, col = 'lightblue')
mtext(side=3, line=1.5, i)
}
for (i in 1:4) {
s2xy <- 25^2 * c(1/i, i)
aligned.pop <- simpopn.bvn(s2xy = s2xy, core = tempgrid, buffer = 100, D = 1, theta = -1)
plotpopn.bvn(aligned.pop, col = 'lightblue')
}
Fig. 1. Elliptical home ranges with varying ratio of major to minor axes as shown. Upper row oriented randomly and independently, lower row with shared random orientation ('randomly aligned').
Simulating elliptical detection data
Generating detection histories
Normally in secr we use sim.capthist to generate capture histories, but that is limited to circular detection functions. The secrBVN function simcapt.bvn is a partial replacement for sim.capthist that models detection with a bivariate normal. Specifically, the cumulative hazard of detection is a constant times the bivariate normal probability density at the detector. The constant is $\lambda_0 2 \pi \sigma_X \sigma_Y$. The user provides a 'popn' object that includes the BVN parameter values ($\sigma_x^2, \sigma_y^2, \theta$) for each animal (row), as generated above by simpopn.bvn. The constant scales the BVN density so that the maximum hazard is $\lambda_0$.
The preceding paragraph describes the default behaviour of simcapt.bvn. If type = uniform is selected then a uniform (flat-topped) elliptical home range is simulated. The uniform probability of detection within the ellipse is controlled by g0, not lambda0.
Wrapper function to generate BVN data and fit SECR model
The function runEllipseSim is a wrapper for the preceding steps (simpopn.bvn, simcapt.bvn) that also fits a standard (circular) SECR model with secr.fit[^footnote1]. The default population has a fixed number of individuals within the rectangular buffered area around the detectors (Ndist = 'fixed').
[^footnote1]: An elliptical SECR model also may be fitted - see Anisotropic home ranges.
For this example we use a $6 \times 6$ grid of binary proximity detectors with 50 metre spacing. The code in secrBVN does not allow for competition among detectors (secr detector type 'multi') or other other secr detector types. A density of 4/ha gives exactly 169 animals in the buffered area. Conditional likelihood is used for speed; the default extractfn = derived is compatible with both CL = TRUE and CL = FALSE. The 200-m buffer allows for the longest ranges ($\sigma_y = 50$ m).
tr <- make.grid(6, 6, spacing = 50, detector = 'proximity')
simrandomBVN36 <- vector('list', 4)
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
simrandomBVN36[[i]] <- runEllipseSim (nrepl = 500, sigmaX, sigmaY, buffer = 200,
ncores = 20, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN',
CL = TRUE, detectfn = 'HHN', details = list(distribution = 'binomial'))
}
The secrBVN function simplot is used to summarize the results
load(file = paste0(simfolder, 'simrandomBVN36.RData'))
par(xpd = FALSE, mar = c(4,4,4,4))
simplot(list(BVN = simrandomBVN36))
Fig. 2. Relative bias of density estimated by fitting circular SECR detection model to elliptical data, with 95\% confidence limit for simulated values.
Next print a summary of the simulation results. There is no apparent effect of range elongation itself on the bias of the estimates.
simsum(list(BVN = simrandomBVN36))
The Appendix has many other examples.
Heterogeneous elliptical home ranges {#heterogeneity}
simpopn.bvn by default generates a population with equal-sized home ranges. The size and elongation of each range may be varied by providing a user-defined function as the argument s2xy:
tr <- make.grid(6,6, spacing = 50, detector = 'proximity')
rs2xy <- function(N, scale = 25) {
aspectratio <- 1 + runif(N) * 3
cbind(scale^2 / aspectratio, scale^2 * aspectratio)
}
pop <- simpopn.bvn(s2xy = rs2xy, core = tr, buffer = 100, D = 1)
par(mfrow = c(1,1), mar = c(4,4,4,4), xpd = TRUE)
plotpopn.bvn(pop, col = 'grey')
Heterogeneous populations may also be simulated in runEllipseSim by passing a function as the argument 'sigmaX'.
sims <- runEllipseSim (10, sigmaX = rs2xy, buffer = 200, ncores = 2,
traps = tr, g0 = 0.2, D = 4, type = 'uniform')
Anisotropic home ranges: a partial solution {#anisotropic}
In principle, we can deal with elongated ranges by replacing Euclidean distances with distances in a space transformed to render home ranges circular (Murphy et al. 2016). The transformation compresses home ranges along their major axis. The transformation uses two parameters: $\psi_A$ for the shared orientation measured in radians and $\psi_R$ for the compression ratio ($\psi_R \ge 1$; $\psi_R = 1$ corresponds to no compression). The parameter $\psi_A$ (psiA) corresponds to the argument theta used by simpopn.bvn. Thanks to Ben Augustine for pointing out the geoR function coords.aniso that lets us do this (Ribeiro and Diggle 2018)[^footnote].
[^footnote]:The package intamap (Pebesma et al. 2010) also offers anisotropic transformation in function rotateAnisotropicData.
The method does not work if the detector array is strictly linear (2-D data are required to estimate the orientation and compression parameters). It is tested in the Appendix on simulated data from a hollow square array.
Distances in transformed space
The code uses a user-defined distance function that computes a Euclidean distance in the transformed space. In this version both transformation parameters are estimated (cf Murphy et al. 2016). Transformation parameters are handled in the undocumented 'miscparm' feature of secr.fit. This allows supernumerary unmodelled parameters to be passed to the distance function. The parameters are included in the vector of coefficients (beta parameters) over which the likelihood is maximised. 'Unmodelled' here means that the parameter takes a single value for all animals, times and places. The link function is 'identity' for $\psi_A$ and 'log' for $\psi_R-1$.
anisodistfn <- function (xy1, xy2, mask) {
if (missing(xy1)) return(character(0))
xy1 <- as.matrix(xy1)
xy2 <- as.matrix(xy2)
miscparm <- attr(mask, 'miscparm')
psiA <- miscparm[1] # anisotropy angle; identity link
psiR <- 1 + exp(miscparm[2]) # anisotropy ratio; log link
aniso.xy1 <- geoR::coords.aniso(xy1, aniso.pars = c(psiA, psiR))
aniso.xy2 <- geoR::coords.aniso(xy2, aniso.pars = c(psiA, psiR))
secr::edist(aniso.xy1, aniso.xy2) # nrow(xy1) x nrow(xy2) matrix
}
Example
Here is a simple application with oblique elliptical home ranges (theta = pi/4).
tr <- make.grid(10, 10, spacing = 25, hollow = TRUE, detector = 'proximity')
pop <- simpopn.bvn(s2xy = (25*c(0.5,2))^2, theta = pi/4, core = tr, buffer = 200,
D = 4, Ndist = 'fixed')
CH <- simcapt.bvn(tr, pop, type = 'BVN', lambda = 0.4, noccasions = 5)
par(mfrow = c(1,1))
plot(CH, tracks = TRUE)
plotpopn.bvn(pop, border='grey', add = TRUE)
plot(tr, add = TRUE)
First we fit a naive circular model.
fit0 <- secr.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE,
details = list(distribution = 'binomial'))
predict(fit0)
Next, the model with transformation to isotropy.
details <- list(distribution = 'binomial', userdist = anisodistfn,
miscparm = c(psiA = 0.5, psiR = 1.5)) # initial values
fit1 <- secr.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE,
details = details)
predict(fit1)
coef(fit1)
The estimated bearing is r round(coef(fit1)['psiA','beta']*360/2/pi,2) degrees. The estimated aspect ratio is r round(1 + exp(coef(fit1)['psiR', 'beta']),2), (95% CI r round(1 + exp(coef(fit1)['psiR', c('lcl','ucl')]), 2)).
Function to fit anisotropic model
The secrBVN function anisotropic.fit streamlines model fitting and does not require anisodistfn to be defined externally. Use it as you would use secr.fit. For example,
fit2 <- anisotropic.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE,
details = list(distribution = "binomial"))
predictAniso(fit2) # estimates of transformation parameters psiA (in degrees) and psiR
Fixed orientation
If the orientation is known then there are two solutions:
- Re-write anisodistfn for miscparm of length 1, with hard-coded value for psiA, or
- Fix the beta parameter corresponding to psiA.
Here we demonstrate (2), fixing psiA at pi/4 (45 degrees). The secr.fit details argument fixedbeta is a vector of values, one per 'beta' parameter, using NA for each beta to be estimated. The first three beta parameters are for D, lambda0 and sigma; the next two are for psiA and psiR. If you model D, lambda0 or sigma there will be extra beta parameters, pushing psiA and psiR back.
fit3 <- anisotropic.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE,
details = list(distribution = "binomial", fixedbeta = c(NA,NA,NA,pi/4,NA)))
predictAniso(fit3) # estimate of psiR only
predict(fit3)
Inadequate data
Detectors in a straight line provide very little information with which to estimate $\psi_A$ and $\psi_R$. Here is an example.
tr <- make.grid(100, 1, spacing = 25, detector = 'proximity')
pop <- simpopn.bvn(s2xy = (25*c(0.5,2))^2, theta = pi/4, core = tr, buffer = 200,
D = 4, Ndist = 'fixed')
CH <- simcapt.bvn(tr, pop, type = 'BVN', lambda = 0.4, noccasions = 5)
fit4 <- anisotropic.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE,
details = list(distribution = "binomial"))
predictAniso(fit4) # estimates of transformation parameters psiA (in degrees) and psiR
predict(fit4)
(Perhaps the parameters are not strictly nonidentifiable).
Package limitations
This package has the limited goal of determining how range elongation affects estimates of density in simple SECR models, and these specific limitations:
- Only binary proximity detectors are supported.
- The spatial distribution of activity centres is assumed to be homogeneous Poisson.
- Ellipses are specified using either 'sigmaX' and 'sigmaY' as separate arguments (
runEllipseSim) or as a vector of the two values, squared ('s2xy' in simpopn.bvn). This is confusing but it's better at this point not to change.
References
Efford, M. G. (2004) Density estimation in live-trapping studies. Oikos 106, 598--610.
Efford, M. G. In prep. Non-circular home ranges and the estimation of population density.
Huggins, R. M. (1989) On the statistical analysis of capture experiments.
Biometrika 76, 133--140.
Ivan, J. S., White, G. C. and Shenk, T. M. (2013) Using simulation to compare methods for
estimating density from capture--recapture data.
Ecology 94, 817--826.
Murphy, S. M., Cox, J. J., Augustine, B. C., Hast, J. T., Guthrie, J. M., Wright, J., McDermott, J., Maehr, S. C. and Plaxico, J. H. (2016) Characterizing recolonization by a reintroduced bear population using genetic spatial capture--recapture. Journal of Wildlife Management 80, 1390--1407.
Pebesma, E., Cornford, D., Dubois, G., Heuvelink, G.B.M., Hristopoulos, D., Pilz, J., Stoehlker, U., Morin, G. and Skoien, J.O. (2010) INTAMAP: the design and implementation of an interoperable automated interpolation web service.
Computers & Geosciences 37, 343--352.
Ribeiro, P. J. Jr and Diggle, P. J. (2018) geoR: Analysis of
Geostatistical Data. R package version 1.7-5.2.1.
https://CRAN.R-project.org/package=geoR.
\pagebreak
Appendix. Code for simulations of Efford (in prep) {#Appendix}
library(secrBVN)
simfolder <- "d:/density communication/noncircularity/paper/simulations/"
nrepl <- 500
ncores <- 20 # for machine with at least 20 cores
runsim <- FALSE # change to TRUE to execute simulations
details <- list(distribution = 'binomial') # all simulations assume fixed N
Functions for heterogeneous aspect ratio.
# Uniform on 1-4
rs2xy <- function(N, scale = 25) {
i <- 1 + runif(N) * 3
cbind(scale^2 /i, scale^2 * i)
}
# 2 classes 1,4
rs2xy2 <- function(N, scale = 25, prob = c(0.5,0.5)) {
i <- sample(c(1,4), size = N, replace = TRUE, prob = prob)
cbind(scale^2 /i, scale^2 * i)
}
Main simulations
10 x 10 grid
Elongated ranges are oriented at random with respect to the grid (default theta = NULL).
tr <- make.grid(10,10, spacing = 50, detector = 'proximity')
simrandom100.3 <- vector('list', 6)
names(simrandom100.3) <- c(1:4,'1-4','1,4')
simrandomBVN100.3 <- simrandom100.3
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, theta = NULL,
D = 4, CL = TRUE, detectfn = 'HHN', details = details, seed = 347)
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
# uniform
args$type <- 'uniform'; args$g0 <- 0.2
simrandom100.3[[i]] <- do.call(runEllipseSim, args)
# bvn
args$type <- 'BVN'; args$lambda0 <- 0.4
simrandomBVN100.3[[i]] <- do.call(runEllipseSim, args)
message ('Completed aspect ratio', i)
}
# heterogeneous aspect ratio 1,4, uniform
args$type <- 'uniform'; args$sigmaX <- rs2xy
simrandom100.3[[5]] <- do.call(runEllipseSim, args)
# heterogeneous aspect ratio 1,4, BVN
args$type <- 'BVN'; args$sigmaX <- rs2xy
simrandomBVN100.3[[5]] <- do.call(runEllipseSim, args)
# heterogeneous aspect ratio 1,4, uniform
args$type <- 'uniform'; args$sigmaX <- rs2xy2
simrandom100.3[[6]] <- do.call(runEllipseSim, args)
#heterogeneous aspect ratio 1,4, BVN
args$type <- 'BVN'; args$sigmaX <- rs2xy2; args$seed <- 347
simrandomBVN100.3[[6]] <- do.call(runEllipseSim, args)
save(simrandom100.3, file = paste0(simfolder, 'simrandom100.3.RData'))
save(simrandomBVN100.3, file = paste0(simfolder, 'simrandomBVN100.3.RData'))
load(file = paste0(simfolder, 'simrandom100.3.RData'))
load(file = paste0(simfolder, 'simrandomBVN100.3.RData'))
par(mfrow = c(1,1), mar = c(4,4,4,4), xpd = FALSE)
# select only first 4 scenarios for plotting
simplot(list(Uniform = simrandom100.3[1:4], BVN = simrandomBVN100.3[1:4]), legend = TRUE)
# tabulate all
simsum(list(Uniform = simrandom100.3, BVN = simrandomBVN100.3), compact = NULL, dec = 4)
simsum(list(Uniform = simrandom100.3, BVN = simrandomBVN100.3), dec = 4)
# compact table for paper
simsum(list(Uniform = simrandom100.3, BVN = simrandomBVN100.3))
# spatial scale parameter
output <- simsum(list(Uniform = simrandom100.3, BVN = simrandomBVN100.3),
component = 'pred', parm = 'sigma', compact = c("av.parmhat", "sd.parmhat"))
lapply(output, '/', 50) # in units of detector spacing
Straight line 36-detectors
tr <- make.grid(36, 1, spacing = 50, detector = 'proximity')
simgridalignlinw1.1 <- vector('list', 6)
names(simgridalignlinw1.1) <- c(1:4,'1-4','1,4')
simgridalignlinw4.1 <- simgridalignlinw3.1 <- simgridalignlinw2.1 <- simgridalignlinw1.1
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN',
CL = TRUE, detectfn = 'HHN', details = details)
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
# bvn perpendicular to traps
args$theta <- 0; simgridalignlinw1.1[[i]] <- do.call(runEllipseSim, args)
# bvn oblique to traps
args$theta <- pi/4; simgridalignlinw2.1[[i]] <- do.call(runEllipseSim, args)
# bvn parallel to traps
args$theta <- pi/2; simgridalignlinw3.1[[i]] <- do.call(runEllipseSim, args)
# bvn random orientation
args$theta <- NULL; simgridalignlinw4.1[[i]] <- do.call(runEllipseSim, args)
}
args <- c(baseargs, list(sigmaX = rs2xy))
args$theta <- 0; simgridalignlinw1.1[[5]] <- do.call(runEllipseSim, args)
args$theta <- pi/4; simgridalignlinw2.1[[5]] <- do.call(runEllipseSim, args)
args$theta <- pi/2; simgridalignlinw3.1[[5]] <- do.call(runEllipseSim, args)
args$theta <- NULL; simgridalignlinw4.1[[5]] <- do.call(runEllipseSim, args)
args <- c(baseargs, list(sigmaX = rs2xy2))
args$theta <- 0; simgridalignlinw1.1[[6]] <- do.call(runEllipseSim, args)
args$theta <- pi/4; simgridalignlinw2.1[[6]] <- do.call(runEllipseSim, args)
args$theta <- pi/2; simgridalignlinw3.1[[6]] <- do.call(runEllipseSim, args)
args$theta <- NULL; simgridalignlinw4.1[[6]] <- do.call(runEllipseSim, args)
save(simgridalignlinw1.1, file = paste0(simfolder, 'simgridalignlinw1.1.RData'))
save(simgridalignlinw2.1, file = paste0(simfolder, 'simgridalignlinw2.1.RData'))
save(simgridalignlinw3.1, file = paste0(simfolder, 'simgridalignlinw3.1.RData'))
save(simgridalignlinw4.1, file = paste0(simfolder, 'simgridalignlinw4.1.RData'))
# compact table linear
load(file = paste0(simfolder, 'simgridalignlinw1.1.RData'))
load(file = paste0(simfolder, 'simgridalignlinw2.1.RData'))
load(file = paste0(simfolder, 'simgridalignlinw3.1.RData'))
load(file = paste0(simfolder, 'simgridalignlinw4.1.RData'))
simsum(list(Perpendicular = simgridalignlinw1.1,
Oblique = simgridalignlinw2.1,
Parallel = simgridalignlinw3.1,
Random = simgridalignlinw4.1))
Hollow square 36 detectors
tr <- make.grid(10, 10, spacing = 50, detector = 'proximity', hollow = TRUE)
simgridalignsqw1.1 <- vector('list', 4) # 'w' for wide
names(simgridalignsqw1.1) <- 1:4
simgridalignsqw3.1 <- simgridalignsqw2.1 <- simgridalignsqw1.1
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN',
CL = TRUE, detectfn = 'HHN', details = details)
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
# bvn aligned to traps
args$theta <- 0; simgridalignsqw1.1[[i]] <- do.call(runEllipseSim, args)
# bvn oblique to traps
args$theta <- pi/4; simgridalignsqw2.1[[i]] <- do.call(runEllipseSim, args)
# bvn random orientation
args$theta <- NULL; simgridalignsqw3.1[[i]] <- do.call(runEllipseSim, args)
}
args <- c(baseargs, list(sigmaX = rs2xy))
args$theta <- 0; simgridalignsqw1.1[[5]] <- do.call(runEllipseSim, args)
args$theta <- pi/4; simgridalignsqw2.1[[5]] <- do.call(runEllipseSim, args)
args$theta <- NULL; simgridalignsqw3.1[[5]] <- do.call(runEllipseSim, args)
args <- c(baseargs, list(sigmaX = rs2xy2))
args$theta <- 0; simgridalignsqw1.1[[6]] <- do.call(runEllipseSim, args)
args$theta <- pi/4; simgridalignsqw2.1[[6]] <- do.call(runEllipseSim, args)
args$theta <- NULL; simgridalignsqw3.1[[6]] <- do.call(runEllipseSim, args)
save(simgridalignsqw1.1, file = paste0(simfolder, 'simgridalignsqw1.1.RData'))
save(simgridalignsqw2.1, file = paste0(simfolder, 'simgridalignsqw2.1.RData'))
save(simgridalignsqw3.1, file = paste0(simfolder, 'simgridalignsqw3.1.RData'))
# compact table square
load(file = paste0(simfolder, 'simgridalignsqw1.1.RData'))
load(file = paste0(simfolder, 'simgridalignsqw2.1.RData'))
load(file = paste0(simfolder, 'simgridalignsqw3.1.RData'))
simsum(list(Aligned = simgridalignsqw1.1,
Oblique = simgridalignsqw2.1,
Random = simgridalignsqw3.1))
Reduced spacing ($\sigma$ instead of $2\sigma$)
10 x 10 grid
tr <- make.grid(10,10, spacing = 25, detector = 'proximity')
simrandomBVNs100.1 <- vector('list', 4)
names(simrandomBVNs100.1) <- 1:4
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, theta = NULL,
D = 4, CL = TRUE, detectfn = 'HHN', details = details, seed = 347,
type = 'BVN', lambda = 0.4)
for (i in 1:4) {
args <- c(baseargs, list(sigmaX = 25/i^0.5, sigmaY = 25*i^0.5))
simrandomBVNs100.1[[i]] <- do.call(runEllipseSim, args)
}
save(simrandomBVNs100.1, file = paste0(simfolder, 'simrandomBVNs100.1.RData'))
Straight line
tr <- make.grid(36, 1, spacing = 25, detector = 'proximity')
simgridalignlinw1.1 <- vector('list', 6)
names(simgridalignlinw1.1) <- c(1:4,'1-4','1,4')
simgridalignlinw4.1 <- simgridalignlinw3.1 <- simgridalignlinw2.1 <- simgridalignlinw1.1
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN', theta = 0,
CL = TRUE, detectfn = 'HHN', details = details)
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
# bvn perpendicular to traps
args$theta <- 0; simgridalignlin1.3[[i]] <- do.call(runEllipseSim, args)
# bvn oblique to traps
args$theta <- pi/4; simgridalignlin2.3[[i]] <- do.call(runEllipseSim, args)
# bvn parallel to traps
args$theta <- pi/2; simgridalignlin3.3[[i]] <- do.call(runEllipseSim, args)
# bvn random orientation
args$theta <- NULL; simgridalignlin4.3[[i]] <- do.call(runEllipseSim, args)
}
args <- c(baseargs, list(sigmaX = rs2xy))
args$theta <- 0; simgridalignlin1.3[[5]] <- do.call(runEllipseSim, args)
args$theta <- pi/4; simgridalignlin2.3[[5]] <- do.call(runEllipseSim, args)
args$theta <- pi/2; simgridalignlin3.3[[5]] <- do.call(runEllipseSim, args)
args$theta <- NULL; simgridalignlin4.3[[5]] <- do.call(runEllipseSim, args)
args <- c(baseargs, list(sigmaX = rs2xy2))
args$theta <- 0; simgridalignlin1.3[[6]] <- do.call(runEllipseSim, args)
args$theta <- pi/4; simgridalignlin2.3[[6]] <- do.call(runEllipseSim, args)
args$theta <- pi/2; simgridalignlin3.3[[6]] <- do.call(runEllipseSim, args)
args$theta <- NULL; simgridalignlin4.3[[6]] <- do.call(runEllipseSim, args)
save(simgridalignlin1.3, file = paste0(simfolder, 'simgridalignlin1.3.RData'))
save(simgridalignlin2.3, file = paste0(simfolder, 'simgridalignlin2.3.RData'))
save(simgridalignlin3.3, file = paste0(simfolder, 'simgridalignlin3.3.RData'))
save(simgridalignlin4.3, file = paste0(simfolder, 'simgridalignlin4.3.RData'))
Hollow square
tr <- make.grid(10, 10, spacing = 25, detector = 'proximity', hollow = TRUE)
simgridalignsq1.1 <- vector('list', 4) # 'w' for wide
names(simgridalignsq1.1) <- 1:4
simgridalignsq3.1 <- simgridalignsq2.1 <- simgridalignsq1.1
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN',
CL = TRUE, detectfn = 'HHN', details = details)
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
# bvn aligned to traps
args$theta <- 0; simgridalignsq1.1[[i]] <- do.call(runEllipseSim, args)
# bvn oblique to traps
args$theta <- pi/4; simgridalignsq2.1[[i]] <- do.call(runEllipseSim, args)
# bvn random orientation
args$theta <- NULL; simgridalignsq3.1[[i]] <- do.call(runEllipseSim, args)
}
args <- c(baseargs, list(sigmaX = rs2xy))
args$theta <- 0; simgridalignsq1.1[[5]] <- do.call(runEllipseSim, args)
args$theta <- pi/4; simgridalignsq2.1[[5]] <- do.call(runEllipseSim, args)
args$theta <- NULL; simgridalignsq3.1[[5]] <- do.call(runEllipseSim, args)
args <- c(baseargs, list(sigmaX = rs2xy2))
args$theta <- 0; simgridalignsq1.1[[6]] <- do.call(runEllipseSim, args)
args$theta <- pi/4; simgridalignsq2.1[[6]] <- do.call(runEllipseSim, args)
args$theta <- NULL; simgridalignsq3.1[[6]] <- do.call(runEllipseSim, args)
save(simgridalignsq1.2, file = paste0(simfolder, 'simgridalignsq1.2.RData'))
save(simgridalignsq2.2, file = paste0(simfolder, 'simgridalignsq2.2.RData'))
save(simgridalignsq3.2, file = paste0(simfolder, 'simgridalignsq3.2.RData'))
Spacing comparisons
10 x 10 grid
load(file = paste0(simfolder, 'simrandomBVN100.3.RData'))
load(file = paste0(simfolder, 'simrandomBVNs100.1.RData'))
par(mfrow=c(1,2))
simplot(list(BVN = simrandomBVNs100.1), legend = FALSE)
mtext(side=3, text = '25 m')
simplot(list(BVN = simrandomBVN100.3[1:4]), legend = FALSE)
mtext(side=3, text = '50 m')
Straight line
load(file = paste0(simfolder, 'simgridalignlin1.3.RData'))
load(file = paste0(simfolder, 'simgridalignlin2.3.RData'))
load(file = paste0(simfolder, 'simgridalignlin3.3.RData'))
load(file = paste0(simfolder, 'simgridalignlin4.3.RData'))
load(file = paste0(simfolder, 'simgridalignlinw1.1.RData'))
load(file = paste0(simfolder, 'simgridalignlinw2.1.RData'))
load(file = paste0(simfolder, 'simgridalignlinw3.1.RData'))
load(file = paste0(simfolder, 'simgridalignlinw4.1.RData'))
par(mfrow=c(1,2))
simplot(list(Perpendicular = simgridalignlin1.3[1:4],
Oblique = simgridalignlin2.3[1:4],
Parallel = simgridalignlin3.3[1:4],
Random = simgridalignlin4.3[1:4]),
legend = TRUE, ylim = c(-1,3.5))
mtext(side=3, text = '25 m')
simplot(list(Perpendicular = simgridalignlinw1.1[1:4],
Oblique = simgridalignlinw2.1[1:4],
Parallel = simgridalignlinw3.1[1:4],
Random = simgridalignlinw4.1[1:4]),
legend = TRUE, ylim = c(-1,3.5))
mtext(side=3, text = '50 m')
Hollow square
load(file = paste0(simfolder, 'simgridalignsq1.2.RData'))
load(file = paste0(simfolder, 'simgridalignsq2.2.RData'))
load(file = paste0(simfolder, 'simgridalignsq3.2.RData'))
load(file = paste0(simfolder, 'simgridalignsqw1.1.RData'))
load(file = paste0(simfolder, 'simgridalignsqw2.1.RData'))
load(file = paste0(simfolder, 'simgridalignsqw3.1.RData'))
par(mfrow =c(1,2))
simplot(list(Aligned = simgridalignsq1.2[1:4],
Oblique = simgridalignsq2.2[1:4],
Random = simgridalignsq3.2[1:4]),
legend = TRUE, ylim = c(-1,3.5))
mtext(side=3, text = '25 m')
simplot(list(Aligned = simgridalignsqw1.1[1:4],
Oblique = simgridalignsqw2.1[1:4],
Random = simgridalignsqw3.1[1:4]),
legend = TRUE, ylim = c(-1,3.5))
mtext(side=3, text = '50 m')
Variations
Low density (0.5/ha)
tr <- make.grid(10, 10, spacing = 50, detector = 'proximity')
simrandom100low.2 <- vector('list', 4)
names(simrandom100low.2) <- 1:4
simrandomBVN100low.2 <- simrandom100low.2
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
D = 0.5, CL = TRUE, detectfn = 'HHN', details = details)
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
args$type <- 'uniform'; args$g0 <- 0.2
simrandom100low.2[[i]] <- do.call(runEllipseSim, args)
args$type <- 'BVN'; args$lambda0 <- 0.4
simrandomBVN100low.2[[i]] <- do.call(runEllipseSim, args)
}
save(simrandom100low.2, file = paste0(simfolder, 'simrandom100low.2.RData'))
save(simrandomBVN100low.2, file = paste0(simfolder, 'simrandomBVN100low.2.RData'))
load(file = paste0(simfolder, 'simrandom100low.2.RData'))
load(file = paste0(simfolder, 'simrandomBVN100low.2.RData'))
par(mfrow = c(1,1), mar = c(4,4,4,4), xpd = FALSE)
simplot(list(Uniform = simrandom100low.2, BVN = simrandomBVN100low.2),
trueval = 0.5, legend = TRUE)
simsum(list(Uniform = simrandom100low.2, BVN = simrandomBVN100low.2),
trueval = 0.5, compact = NULL)
Small array (6 x 6 grid)
tr <- make.grid(6, 6, spacing = 50, detector = 'proximity')
simrandom36.2 <- vector('list', 4)
names(simrandom36.2) <- 1:4
simrandomBVN36.2 <- simrandom36.2
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
D = 4, CL = TRUE, detectfn = 'HHN', details = details)
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
args$type <- 'uniform'; args$g0 <- 0.2
simrandom36.2[[i]] <- do.call(runEllipseSim, args)
args$type <- 'BVN'; args$lambda0 <- 0.4
simrandomBVN36.2[[i]] <- do.call(runEllipseSim, args)
}
save(simrandom36.2, file = paste0(simfolder, 'simrandom36.2.RData'))
save(simrandomBVN36.2, file = paste0(simfolder, 'simrandomBVN36.2.RData'))
load(file = paste0(simfolder, 'simrandom36.2.RData'))
load(file = paste0(simfolder, 'simrandomBVN36.2.RData'))
par(mfrow = c(1,1), mar = c(4,4,4,4), xpd = FALSE)
simplot(list(Uniform = simrandom36.2[1:4], BVN = simrandomBVN36.2[1:4]), legend = TRUE)
simsum(list(Uniform = simrandom36.2, BVN = simrandomBVN36.2))
Common random orientation
tr <- make.grid(10, 10, spacing = 50, detector = 'proximity')
simrandom100C.2 <- vector('list', 4)
names(simrandom100C.2) <- 1:4
simrandomBVN100C.2 <- simrandom100C.2
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
theta = -1, D = 4, CL = TRUE, detectfn = 'HHN', details = details)
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
args$type <- 'uniform'; args$g0 <- 0.2
simrandom100C.2[[i]] <- do.call(runEllipseSim, args)
args$type <- 'BVN'; args$lambda0 <- 0.4
simrandomBVN100C.2[[i]] <- do.call(runEllipseSim, args)
}
save(simrandom100C.2, file = paste0(simfolder, 'simrandom100C.2.RData'))
save(simrandomBVN100C.2, file = paste0(simfolder, 'simrandomBVN100C.2.RData'))
load(file = paste0(simfolder, 'simrandom100C.2.RData'))
load(file = paste0(simfolder, 'simrandomBVN100C.2.RData'))
simplot(list(Uniform = simrandom100C.2, BVN = simrandomBVN100C.2), legend = TRUE)
Grid with varying orientation
Spacings narrow and wide (w); orientation parallel (p) and oblique(o).
Only consider aspect ratio 4 for now.
tr <- make.grid(10,10, spacing = 25, detector = 'proximity')
simrandomBVN100p.1 <- vector('list', 4)
names(simrandomBVN100p.1) <- 1:4
simrandomBVN100o.1 <- simrandomBVN100p.1
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
type = 'BVN', lambda0 = 0.4,
D = 4, CL = TRUE, detectfn = 'HHN', details = details)
for (i in 4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY, theta = 0))
simrandomBVN100p.1[[i]] <- do.call(runEllipseSim, args)
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY, theta = pi/4))
simrandomBVN100o.1[[i]] <- do.call(runEllipseSim, args)
}
tr <- make.grid(10,10, spacing = 50, detector = 'proximity')
simrandomBVN100pw.1 <- vector('list', 4)
names(simrandomBVN100pw.1) <- 1:4
simrandomBVN100ow.1 <- simrandomBVN100pw.1
for (i in 4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY, theta = 0))
simrandomBVN100pw.1[[i]] <- do.call(runEllipseSim, args)
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY, theta = pi/4))
simrandomBVN100ow.1[[i]] <- do.call(runEllipseSim, args)
}
save(simrandomBVN100p.1, file = paste0(simfolder, 'simrandomBVN100p.1.RData'))
save(simrandomBVN100o.1, file = paste0(simfolder, 'simrandomBVN100o.1.RData'))
save(simrandomBVN100pw.1, file = paste0(simfolder, 'simrandomBVN100pw.1.RData'))
save(simrandomBVN100ow.1, file = paste0(simfolder, 'simrandomBVN100ow.1.RData'))
Anisotropic model
Using function anisotropic.fit. Extract the fitted coefficients corresponding to psiA and psiR with the coef method for secr objects, and rely on direct estimation of density (CL = FALSE).
Hollow square
tr <- make.grid(10, 10, spacing = 25, detector = 'proximity', hollow = TRUE)
simaniso1 <- vector('list', 4); names(simaniso1) <- 1:4
simaniso4 <- simaniso3 <- simaniso2 <- simaniso1
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN',
details = list(distribution = 'binomial'), extractfn = coef,
secrfn = 'anisotropic.fit')
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
args$theta <- 0 # aligned to detectors
simaniso1[[i]] <- do.call(runEllipseSim, args)
args$theta <- pi/4 # oblique to detectors
simaniso2[[i]] <- do.call(runEllipseSim, args)
args$theta <- NULL # random orientation
simaniso3[[i]] <- do.call(runEllipseSim, args)
args$theta <- -1 # common random orientation
simaniso4[[i]] <- do.call(runEllipseSim, args)
}
save(simaniso1, file = paste0(simfolder, 'simaniso1.RData'))
save(simaniso2, file = paste0(simfolder, 'simaniso2.RData'))
save(simaniso3, file = paste0(simfolder, 'simaniso3.RData'))
save(simaniso4, file = paste0(simfolder, 'simaniso4.RData'))
load(file = paste0(simfolder, 'simgridalignsq1.2.RData'))
load(file = paste0(simfolder, 'simgridalignsq2.2.RData'))
load(file = paste0(simfolder, 'simgridalignsq3.2.RData'))
load(file = paste0(simfolder, 'simaniso1.RData'))
load(file = paste0(simfolder, 'simaniso2.RData'))
load(file = paste0(simfolder, 'simaniso3.RData'))
load(file = paste0(simfolder, 'simaniso4.RData'))
par(mfrow = c(1,2))
simplot(list(Aligned = simgridalignsq1.2[1:4],
Oblique = simgridalignsq2.2[1:4],
Random = simgridalignsq3.2[1:4]),
ylim = c(-0.4,0.45), legend = TRUE, pchi = c(23, 22, 16))
text(-0.4, 0.5, 'a.', cex = 1.45, xpd = TRUE)
text(3, -0.35, 'isotropic')
simplot(list(Aligned = simaniso1,
Oblique = simaniso2,
Random = simaniso3),
component = "pred",
ylim = c(-0.4,0.45), legend = TRUE, pchi = c(23, 22, 16))
text(-0.4, 0.5, 'b.', cex = 1.45, xpd = TRUE)
text(3, -0.35, 'anisotropic')
Tabular summary for anisotropic model.
simsum(list(Aligned = simaniso1,
Oblique = simaniso2,
Random = simaniso3),
component = "pred")
Tabular comparison of isotropic and anisotropic models for hollow grid.
tab1 <- simsum(list(Aligned = simgridalignsq1.2[1:4],
Oblique = simgridalignsq2.2[1:4],
Random = simgridalignsq3.2[1:4])) # find D-hat in the 'fit' component
tab2 <- simsum(list(Aligned = simaniso1,
Oblique = simaniso2,
Random = simaniso3),
component = "pred") # find D-hat in the 'pred' component
fn <- function(t1,t2) cbind(t1[,-1],t2[,-1])
mapply(fn, tab1, tab2, SIMPLIFY = FALSE)
Straight line
tr <- make.grid(36, 1, spacing = 25, detector = 'proximity')
simanisolin1 <- vector('list', 4); names(simanisolin1) <- 1:4
simanisolin5 <- simanisolin4 <- simanisolin3 <- simanisolin2 <- simanisolin1
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN',
details = list(distribution = 'binomial'), extractfn = coef,
secrfn = 'anisotropic.fit')
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
args$theta <- 0 # perpendicular to detectors
simanisolin1[[i]] <- do.call(runEllipseSim, args)
args$theta <- pi/4 # oblique to detectors
simanisolin2[[i]] <- do.call(runEllipseSim, args)
args$theta <- NULL # random orientation
simanisolin3[[i]] <- do.call(runEllipseSim, args)
args$theta <- -1 # common random orientation
simanisolin4[[i]] <- do.call(runEllipseSim, args)
args$theta <- pi/2 # parallel
simanisolin5[[i]] <- do.call(runEllipseSim, args) # partial run 2018-08-23
}
save(simanisolin1, file = paste0(simfolder, 'simanisolin1.RData'))
save(simanisolin2, file = paste0(simfolder, 'simanisolin2.RData'))
save(simanisolin3, file = paste0(simfolder, 'simanisolin3.RData'))
save(simanisolin4, file = paste0(simfolder, 'simanisolin4.RData'))
save(simanisolin5, file = paste0(simfolder, 'simanisolin5.RData'))
load(file = paste0(simfolder, 'simgridalignlin1.3.RData'))
load(file = paste0(simfolder, 'simgridalignlin2.3.RData'))
load(file = paste0(simfolder, 'simgridalignlin4.3.RData'))
load(file = paste0(simfolder, 'simanisolin1.RData'))
load(file = paste0(simfolder, 'simanisolin2.RData'))
load(file = paste0(simfolder, 'simanisolin3.RData'))
par(mfrow = c(1,2))
simplot(list(Aligned = simgridalignlin1.3[1:4],
Oblique = simgridalignlin2.3[1:4],
Random = simgridalignlin4.3[1:4]),
ylim = c(-1,3.5), legend = TRUE, pchi = c(23, 22, 16))
text(-0.1, 3.5, 'a.', cex = 1.45, xpd = TRUE)
text(1, -0.8, 'isotropic', adj=0)
simplot(list(Perpendicular = simanisolin1,
Oblique = simanisolin2,
Random = simanisolin3),
component = "pred",
ylim = c(-1,3.5), legend = TRUE, pchi = c(23, 22, 16))
text(-0.1, 3.5, 'b.', cex = 1.45, xpd = TRUE)
text(1, -0.8, 'anisotropic', adj=0)
Conclude that anisotropic transformation largely ineffective with linear array.
10 x 10 grid
tr <- make.grid(10, 10, spacing = 25, detector = 'proximity')
simanisogrid1 <- vector('list', 4); names(simanisogrid1) <- 1:4
simanisogrid4 <- simanisogrid3 <- simanisogrid2 <- simanisogrid1
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN',
details = list(distribution = 'binomial'), extractfn = coef,
secrfn = 'anisotropic.fit')
for (i in 1:4) {
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY))
args$theta <- 0 # aligned to detectors
simanisogrid1[[i]] <- do.call(runEllipseSim, args)
args$theta <- pi/4 # oblique to detectors
simanisogrid2[[i]] <- do.call(runEllipseSim, args)
args$theta <- NULL # random orientation
simanisogrid3[[i]] <- do.call(runEllipseSim, args)
args$theta <- -1 # common random orientation
simanisogrid4[[i]] <- do.call(runEllipseSim, args)
}
save(simanisogrid1, file = paste0(simfolder, 'simanisogrid1.RData'))
save(simanisogrid2, file = paste0(simfolder, 'simanisogrid2.RData'))
save(simanisogrid3, file = paste0(simfolder, 'simanisogrid3.RData'))
save(simanisogrid4, file = paste0(simfolder, 'simanisogrid4.RData'))
load(file = paste0(simfolder, 'simanisogrid3.RData'))
load(file = paste0(simfolder, 'simrandomBVNs100.1.RData'))
simsum(list(Simple = simrandomBVNs100.1))
simsum(list(Anisotropic = simanisogrid3), component = "pred")
Performs well with full grid data (narrow spacing).
Now run all geometries with wide spacing (w) for aspect ratio 4.0.
tr <- make.grid(10, 10, spacing = 50, detector = 'proximity')
simanisogridw1 <- vector('list', 4); names(simanisogridw1) <- 1:4
simanisogridw3 <- simanisogridw2 <- simanisogridw1
i <- 4
sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN',
details = list(distribution = 'binomial'), extractfn = coef,
secrfn = 'anisotropic.fit', sigmaX = sigmaX, sigmaY = sigmaY)
args <- c(baseargs, list(theta = 0)) # aligned to detectors
simanisogridw1[[i]] <- do.call(runEllipseSim, args)
args$theta <- pi/4 # oblique to detectors
simanisogridw2[[i]] <- do.call(runEllipseSim, args)
args$theta <- NULL # random orientation
simanisogridw3[[i]] <- do.call(runEllipseSim, args)
save(simanisogridw1, file = paste0(simfolder, 'simanisogridw1.RData'))
save(simanisogridw2, file = paste0(simfolder, 'simanisogridw2.RData'))
save(simanisogridw3, file = paste0(simfolder, 'simanisogridw3.RData'))
tr <- make.grid(36, 1, spacing = 50, detector = 'proximity')
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN',
details = list(distribution = 'binomial'), extractfn = coef,
secrfn = 'anisotropic.fit', sigmaX = sigmaX, sigmaY = sigmaY)
simanisolinw1 <- vector('list', 4); names(simanisolinw1) <- 1:4
simanisolinw4 <- simanisolinw3 <- simanisolinw2 <- simanisolinw1
args <- c(baseargs, list(theta = 0)) # perpendicular to detectors
simanisolinw1[[i]] <- do.call(runEllipseSim, args)
args$theta <- pi/4 # oblique to detectors
simanisolinw2[[i]] <- do.call(runEllipseSim, args)
args$theta <- pi/2 # parallel to detectors
simanisolinw3[[i]] <- do.call(runEllipseSim, args)
args$theta <- NULL # random orientation
simanisolinw4[[i]] <- do.call(runEllipseSim, args)
save(simanisolinw1, file = paste0(simfolder, 'simanisolinw1.RData'))
save(simanisolinw2, file = paste0(simfolder, 'simanisolinw2.RData'))
save(simanisolinw3, file = paste0(simfolder, 'simanisolinw3.RData'))
save(simanisolinw4, file = paste0(simfolder, 'simanisolinw4.RData'))
tr <- make.grid(10, 10, spacing = 50, detector = 'proximity', hollow = TRUE)
simanisow1 <- vector('list', 4); names(simanisow1) <- 1:4
simanisow3 <- simanisow2 <- simanisow1
baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr,
lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN',
details = list(distribution = 'binomial'), extractfn = coef,
secrfn = 'anisotropic.fit', sigmaX = sigmaX, sigmaY = sigmaY)
args <- c(baseargs, list(theta = 0)) # aligned to detectors
simanisow1[[i]] <- do.call(runEllipseSim, args)
args$theta <- pi/4 # oblique to detectors
simanisow2[[i]] <- do.call(runEllipseSim, args)
args$theta <- NULL # random orientation
simanisow3[[i]] <- do.call(runEllipseSim, args)
save(simanisow1, file = paste0(simfolder, 'simanisow1.RData'))
save(simanisow2, file = paste0(simfolder, 'simanisow2.RData'))
save(simanisow3, file = paste0(simfolder, 'simanisow3.RData'))
MurrayEfford/secrBVN documentation built on Oct. 14, 2022, 5:20 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
\vspace{16pt}
The small package secrBVN is used to evaluate the performance of SECR estimators when home ranges are bivariate normal (BVN) or uniform (flat-topped) ellipses. We assume detection hazard is directly proportional to home range utilisation (activity). Code to use secrBVN for the simulations of Efford (in prep.) is provided in the Appendix.
The key user-visible functions are simpopn.bvn, plotpopn.bvn, simcapt.bvn, runEllipseSim, simsum, simplot and anisotropic.fit.
\vspace{12pt}
Generating and plotting elliptical home ranges
simpopn.bvn is a wrapper for the secr function sim.popn that adds attributes specifying a bivariate normal home range shape, size and orientation for each individual. By default, shape and size are the same for all individuals, but a mechanism is provided to vary them individually (see Heterogeneous elliptical home ranges).
First, load the package.
library(secrBVN) simfolder <- "d:/density communication/noncircularity/paper/simulations/" runall <- FALSE # skip lengthy simulations
tempgrid <- make.grid(nx = 10, ny = 10) par(mfrow = c(2,4), mar = c(2,2,2.6,2), xpd = TRUE) for (i in 1:4) { s2xy <- 25^2 * c(1/i, i) random.pop <- simpopn.bvn(s2xy = s2xy, core = tempgrid, buffer = 100, D = 1) plotpopn.bvn(random.pop, col = 'lightblue') mtext(side=3, line=1.5, i) } for (i in 1:4) { s2xy <- 25^2 * c(1/i, i) aligned.pop <- simpopn.bvn(s2xy = s2xy, core = tempgrid, buffer = 100, D = 1, theta = -1) plotpopn.bvn(aligned.pop, col = 'lightblue') }
Fig. 1. Elliptical home ranges with varying ratio of major to minor axes as shown. Upper row oriented randomly and independently, lower row with shared random orientation ('randomly aligned').
Simulating elliptical detection data
Generating detection histories
Normally in secr we use sim.capthist to generate capture histories, but that is limited to circular detection functions. The secrBVN function simcapt.bvn is a partial replacement for sim.capthist that models detection with a bivariate normal. Specifically, the cumulative hazard of detection is a constant times the bivariate normal probability density at the detector. The constant is $\lambda_0 2 \pi \sigma_X \sigma_Y$. The user provides a 'popn' object that includes the BVN parameter values ($\sigma_x^2, \sigma_y^2, \theta$) for each animal (row), as generated above by simpopn.bvn. The constant scales the BVN density so that the maximum hazard is $\lambda_0$.
The preceding paragraph describes the default behaviour of simcapt.bvn. If type = uniform is selected then a uniform (flat-topped) elliptical home range is simulated. The uniform probability of detection within the ellipse is controlled by g0, not lambda0.
Wrapper function to generate BVN data and fit SECR model
The function runEllipseSim is a wrapper for the preceding steps (simpopn.bvn, simcapt.bvn) that also fits a standard (circular) SECR model with secr.fit[^footnote1]. The default population has a fixed number of individuals within the rectangular buffered area around the detectors (Ndist = 'fixed').
[^footnote1]: An elliptical SECR model also may be fitted - see Anisotropic home ranges.
For this example we use a $6 \times 6$ grid of binary proximity detectors with 50 metre spacing. The code in secrBVN does not allow for competition among detectors (secr detector type 'multi') or other other secr detector types. A density of 4/ha gives exactly 169 animals in the buffered area. Conditional likelihood is used for speed; the default extractfn = derived is compatible with both CL = TRUE and CL = FALSE. The 200-m buffer allows for the longest ranges ($\sigma_y = 50$ m).
tr <- make.grid(6, 6, spacing = 50, detector = 'proximity') simrandomBVN36 <- vector('list', 4) for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 simrandomBVN36[[i]] <- runEllipseSim (nrepl = 500, sigmaX, sigmaY, buffer = 200, ncores = 20, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = TRUE, detectfn = 'HHN', details = list(distribution = 'binomial')) }
The secrBVN function simplot is used to summarize the results
load(file = paste0(simfolder, 'simrandomBVN36.RData'))
par(xpd = FALSE, mar = c(4,4,4,4)) simplot(list(BVN = simrandomBVN36))
Fig. 2. Relative bias of density estimated by fitting circular SECR detection model to elliptical data, with 95\% confidence limit for simulated values.
Next print a summary of the simulation results. There is no apparent effect of range elongation itself on the bias of the estimates.
simsum(list(BVN = simrandomBVN36))
The Appendix has many other examples.
Heterogeneous elliptical home ranges {#heterogeneity}
simpopn.bvn by default generates a population with equal-sized home ranges. The size and elongation of each range may be varied by providing a user-defined function as the argument s2xy:
tr <- make.grid(6,6, spacing = 50, detector = 'proximity') rs2xy <- function(N, scale = 25) { aspectratio <- 1 + runif(N) * 3 cbind(scale^2 / aspectratio, scale^2 * aspectratio) } pop <- simpopn.bvn(s2xy = rs2xy, core = tr, buffer = 100, D = 1) par(mfrow = c(1,1), mar = c(4,4,4,4), xpd = TRUE) plotpopn.bvn(pop, col = 'grey')
Heterogeneous populations may also be simulated in runEllipseSim by passing a function as the argument 'sigmaX'.
sims <- runEllipseSim (10, sigmaX = rs2xy, buffer = 200, ncores = 2, traps = tr, g0 = 0.2, D = 4, type = 'uniform')
Anisotropic home ranges: a partial solution {#anisotropic}
In principle, we can deal with elongated ranges by replacing Euclidean distances with distances in a space transformed to render home ranges circular (Murphy et al. 2016). The transformation compresses home ranges along their major axis. The transformation uses two parameters: $\psi_A$ for the shared orientation measured in radians and $\psi_R$ for the compression ratio ($\psi_R \ge 1$; $\psi_R = 1$ corresponds to no compression). The parameter $\psi_A$ (psiA) corresponds to the argument theta used by simpopn.bvn. Thanks to Ben Augustine for pointing out the geoR function coords.aniso that lets us do this (Ribeiro and Diggle 2018)[^footnote].
[^footnote]:The package intamap (Pebesma et al. 2010) also offers anisotropic transformation in function rotateAnisotropicData.
The method does not work if the detector array is strictly linear (2-D data are required to estimate the orientation and compression parameters). It is tested in the Appendix on simulated data from a hollow square array.
Distances in transformed space
The code uses a user-defined distance function that computes a Euclidean distance in the transformed space. In this version both transformation parameters are estimated (cf Murphy et al. 2016). Transformation parameters are handled in the undocumented 'miscparm' feature of secr.fit. This allows supernumerary unmodelled parameters to be passed to the distance function. The parameters are included in the vector of coefficients (beta parameters) over which the likelihood is maximised. 'Unmodelled' here means that the parameter takes a single value for all animals, times and places. The link function is 'identity' for $\psi_A$ and 'log' for $\psi_R-1$.
anisodistfn <- function (xy1, xy2, mask) { if (missing(xy1)) return(character(0)) xy1 <- as.matrix(xy1) xy2 <- as.matrix(xy2) miscparm <- attr(mask, 'miscparm') psiA <- miscparm[1] # anisotropy angle; identity link psiR <- 1 + exp(miscparm[2]) # anisotropy ratio; log link aniso.xy1 <- geoR::coords.aniso(xy1, aniso.pars = c(psiA, psiR)) aniso.xy2 <- geoR::coords.aniso(xy2, aniso.pars = c(psiA, psiR)) secr::edist(aniso.xy1, aniso.xy2) # nrow(xy1) x nrow(xy2) matrix }
Example
Here is a simple application with oblique elliptical home ranges (theta = pi/4).
tr <- make.grid(10, 10, spacing = 25, hollow = TRUE, detector = 'proximity') pop <- simpopn.bvn(s2xy = (25*c(0.5,2))^2, theta = pi/4, core = tr, buffer = 200, D = 4, Ndist = 'fixed') CH <- simcapt.bvn(tr, pop, type = 'BVN', lambda = 0.4, noccasions = 5)
par(mfrow = c(1,1)) plot(CH, tracks = TRUE) plotpopn.bvn(pop, border='grey', add = TRUE) plot(tr, add = TRUE)
First we fit a naive circular model.
fit0 <- secr.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE, details = list(distribution = 'binomial')) predict(fit0)
Next, the model with transformation to isotropy.
details <- list(distribution = 'binomial', userdist = anisodistfn, miscparm = c(psiA = 0.5, psiR = 1.5)) # initial values fit1 <- secr.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE, details = details) predict(fit1) coef(fit1)
The estimated bearing is r round(coef(fit1)['psiA','beta']*360/2/pi,2) degrees. The estimated aspect ratio is r round(1 + exp(coef(fit1)['psiR', 'beta']),2), (95% CI r round(1 + exp(coef(fit1)['psiR', c('lcl','ucl')]), 2)).
Function to fit anisotropic model
The secrBVN function anisotropic.fit streamlines model fitting and does not require anisodistfn to be defined externally. Use it as you would use secr.fit. For example,
fit2 <- anisotropic.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE, details = list(distribution = "binomial")) predictAniso(fit2) # estimates of transformation parameters psiA (in degrees) and psiR
Fixed orientation
If the orientation is known then there are two solutions:
- Re-write anisodistfn for miscparm of length 1, with hard-coded value for psiA, or
- Fix the beta parameter corresponding to psiA.
Here we demonstrate (2), fixing psiA at pi/4 (45 degrees). The secr.fit details argument fixedbeta is a vector of values, one per 'beta' parameter, using NA for each beta to be estimated. The first three beta parameters are for D, lambda0 and sigma; the next two are for psiA and psiR. If you model D, lambda0 or sigma there will be extra beta parameters, pushing psiA and psiR back.
fit3 <- anisotropic.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE, details = list(distribution = "binomial", fixedbeta = c(NA,NA,NA,pi/4,NA))) predictAniso(fit3) # estimate of psiR only predict(fit3)
Inadequate data
Detectors in a straight line provide very little information with which to estimate $\psi_A$ and $\psi_R$. Here is an example.
tr <- make.grid(100, 1, spacing = 25, detector = 'proximity') pop <- simpopn.bvn(s2xy = (25*c(0.5,2))^2, theta = pi/4, core = tr, buffer = 200, D = 4, Ndist = 'fixed') CH <- simcapt.bvn(tr, pop, type = 'BVN', lambda = 0.4, noccasions = 5) fit4 <- anisotropic.fit(CH, buffer = 200, detectfn = 'HHN', trace = FALSE, details = list(distribution = "binomial")) predictAniso(fit4) # estimates of transformation parameters psiA (in degrees) and psiR predict(fit4)
(Perhaps the parameters are not strictly nonidentifiable).
Package limitations
This package has the limited goal of determining how range elongation affects estimates of density in simple SECR models, and these specific limitations:
- Only binary proximity detectors are supported.
- The spatial distribution of activity centres is assumed to be homogeneous Poisson.
- Ellipses are specified using either 'sigmaX' and 'sigmaY' as separate arguments (
runEllipseSim) or as a vector of the two values, squared ('s2xy' insimpopn.bvn). This is confusing but it's better at this point not to change.
References
Efford, M. G. (2004) Density estimation in live-trapping studies. Oikos 106, 598--610.
Efford, M. G. In prep. Non-circular home ranges and the estimation of population density.
Huggins, R. M. (1989) On the statistical analysis of capture experiments. Biometrika 76, 133--140.
Ivan, J. S., White, G. C. and Shenk, T. M. (2013) Using simulation to compare methods for estimating density from capture--recapture data. Ecology 94, 817--826.
Murphy, S. M., Cox, J. J., Augustine, B. C., Hast, J. T., Guthrie, J. M., Wright, J., McDermott, J., Maehr, S. C. and Plaxico, J. H. (2016) Characterizing recolonization by a reintroduced bear population using genetic spatial capture--recapture. Journal of Wildlife Management 80, 1390--1407.
Pebesma, E., Cornford, D., Dubois, G., Heuvelink, G.B.M., Hristopoulos, D., Pilz, J., Stoehlker, U., Morin, G. and Skoien, J.O. (2010) INTAMAP: the design and implementation of an interoperable automated interpolation web service. Computers & Geosciences 37, 343--352.
Ribeiro, P. J. Jr and Diggle, P. J. (2018) geoR: Analysis of Geostatistical Data. R package version 1.7-5.2.1. https://CRAN.R-project.org/package=geoR.
\pagebreak
Appendix. Code for simulations of Efford (in prep) {#Appendix}
library(secrBVN) simfolder <- "d:/density communication/noncircularity/paper/simulations/" nrepl <- 500 ncores <- 20 # for machine with at least 20 cores runsim <- FALSE # change to TRUE to execute simulations details <- list(distribution = 'binomial') # all simulations assume fixed N
Functions for heterogeneous aspect ratio.
# Uniform on 1-4 rs2xy <- function(N, scale = 25) { i <- 1 + runif(N) * 3 cbind(scale^2 /i, scale^2 * i) } # 2 classes 1,4 rs2xy2 <- function(N, scale = 25, prob = c(0.5,0.5)) { i <- sample(c(1,4), size = N, replace = TRUE, prob = prob) cbind(scale^2 /i, scale^2 * i) }
Main simulations
10 x 10 grid
Elongated ranges are oriented at random with respect to the grid (default theta = NULL).
tr <- make.grid(10,10, spacing = 50, detector = 'proximity') simrandom100.3 <- vector('list', 6) names(simrandom100.3) <- c(1:4,'1-4','1,4') simrandomBVN100.3 <- simrandom100.3 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, theta = NULL, D = 4, CL = TRUE, detectfn = 'HHN', details = details, seed = 347) for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) # uniform args$type <- 'uniform'; args$g0 <- 0.2 simrandom100.3[[i]] <- do.call(runEllipseSim, args) # bvn args$type <- 'BVN'; args$lambda0 <- 0.4 simrandomBVN100.3[[i]] <- do.call(runEllipseSim, args) message ('Completed aspect ratio', i) } # heterogeneous aspect ratio 1,4, uniform args$type <- 'uniform'; args$sigmaX <- rs2xy simrandom100.3[[5]] <- do.call(runEllipseSim, args) # heterogeneous aspect ratio 1,4, BVN args$type <- 'BVN'; args$sigmaX <- rs2xy simrandomBVN100.3[[5]] <- do.call(runEllipseSim, args) # heterogeneous aspect ratio 1,4, uniform args$type <- 'uniform'; args$sigmaX <- rs2xy2 simrandom100.3[[6]] <- do.call(runEllipseSim, args) #heterogeneous aspect ratio 1,4, BVN args$type <- 'BVN'; args$sigmaX <- rs2xy2; args$seed <- 347 simrandomBVN100.3[[6]] <- do.call(runEllipseSim, args) save(simrandom100.3, file = paste0(simfolder, 'simrandom100.3.RData')) save(simrandomBVN100.3, file = paste0(simfolder, 'simrandomBVN100.3.RData'))
load(file = paste0(simfolder, 'simrandom100.3.RData')) load(file = paste0(simfolder, 'simrandomBVN100.3.RData')) par(mfrow = c(1,1), mar = c(4,4,4,4), xpd = FALSE) # select only first 4 scenarios for plotting simplot(list(Uniform = simrandom100.3[1:4], BVN = simrandomBVN100.3[1:4]), legend = TRUE) # tabulate all simsum(list(Uniform = simrandom100.3, BVN = simrandomBVN100.3), compact = NULL, dec = 4) simsum(list(Uniform = simrandom100.3, BVN = simrandomBVN100.3), dec = 4)
# compact table for paper simsum(list(Uniform = simrandom100.3, BVN = simrandomBVN100.3)) # spatial scale parameter output <- simsum(list(Uniform = simrandom100.3, BVN = simrandomBVN100.3), component = 'pred', parm = 'sigma', compact = c("av.parmhat", "sd.parmhat")) lapply(output, '/', 50) # in units of detector spacing
Straight line 36-detectors
tr <- make.grid(36, 1, spacing = 50, detector = 'proximity') simgridalignlinw1.1 <- vector('list', 6) names(simgridalignlinw1.1) <- c(1:4,'1-4','1,4') simgridalignlinw4.1 <- simgridalignlinw3.1 <- simgridalignlinw2.1 <- simgridalignlinw1.1 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = TRUE, detectfn = 'HHN', details = details) for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) # bvn perpendicular to traps args$theta <- 0; simgridalignlinw1.1[[i]] <- do.call(runEllipseSim, args) # bvn oblique to traps args$theta <- pi/4; simgridalignlinw2.1[[i]] <- do.call(runEllipseSim, args) # bvn parallel to traps args$theta <- pi/2; simgridalignlinw3.1[[i]] <- do.call(runEllipseSim, args) # bvn random orientation args$theta <- NULL; simgridalignlinw4.1[[i]] <- do.call(runEllipseSim, args) } args <- c(baseargs, list(sigmaX = rs2xy)) args$theta <- 0; simgridalignlinw1.1[[5]] <- do.call(runEllipseSim, args) args$theta <- pi/4; simgridalignlinw2.1[[5]] <- do.call(runEllipseSim, args) args$theta <- pi/2; simgridalignlinw3.1[[5]] <- do.call(runEllipseSim, args) args$theta <- NULL; simgridalignlinw4.1[[5]] <- do.call(runEllipseSim, args) args <- c(baseargs, list(sigmaX = rs2xy2)) args$theta <- 0; simgridalignlinw1.1[[6]] <- do.call(runEllipseSim, args) args$theta <- pi/4; simgridalignlinw2.1[[6]] <- do.call(runEllipseSim, args) args$theta <- pi/2; simgridalignlinw3.1[[6]] <- do.call(runEllipseSim, args) args$theta <- NULL; simgridalignlinw4.1[[6]] <- do.call(runEllipseSim, args) save(simgridalignlinw1.1, file = paste0(simfolder, 'simgridalignlinw1.1.RData')) save(simgridalignlinw2.1, file = paste0(simfolder, 'simgridalignlinw2.1.RData')) save(simgridalignlinw3.1, file = paste0(simfolder, 'simgridalignlinw3.1.RData')) save(simgridalignlinw4.1, file = paste0(simfolder, 'simgridalignlinw4.1.RData'))
# compact table linear load(file = paste0(simfolder, 'simgridalignlinw1.1.RData')) load(file = paste0(simfolder, 'simgridalignlinw2.1.RData')) load(file = paste0(simfolder, 'simgridalignlinw3.1.RData')) load(file = paste0(simfolder, 'simgridalignlinw4.1.RData')) simsum(list(Perpendicular = simgridalignlinw1.1, Oblique = simgridalignlinw2.1, Parallel = simgridalignlinw3.1, Random = simgridalignlinw4.1))
Hollow square 36 detectors
tr <- make.grid(10, 10, spacing = 50, detector = 'proximity', hollow = TRUE) simgridalignsqw1.1 <- vector('list', 4) # 'w' for wide names(simgridalignsqw1.1) <- 1:4 simgridalignsqw3.1 <- simgridalignsqw2.1 <- simgridalignsqw1.1 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = TRUE, detectfn = 'HHN', details = details) for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) # bvn aligned to traps args$theta <- 0; simgridalignsqw1.1[[i]] <- do.call(runEllipseSim, args) # bvn oblique to traps args$theta <- pi/4; simgridalignsqw2.1[[i]] <- do.call(runEllipseSim, args) # bvn random orientation args$theta <- NULL; simgridalignsqw3.1[[i]] <- do.call(runEllipseSim, args) } args <- c(baseargs, list(sigmaX = rs2xy)) args$theta <- 0; simgridalignsqw1.1[[5]] <- do.call(runEllipseSim, args) args$theta <- pi/4; simgridalignsqw2.1[[5]] <- do.call(runEllipseSim, args) args$theta <- NULL; simgridalignsqw3.1[[5]] <- do.call(runEllipseSim, args) args <- c(baseargs, list(sigmaX = rs2xy2)) args$theta <- 0; simgridalignsqw1.1[[6]] <- do.call(runEllipseSim, args) args$theta <- pi/4; simgridalignsqw2.1[[6]] <- do.call(runEllipseSim, args) args$theta <- NULL; simgridalignsqw3.1[[6]] <- do.call(runEllipseSim, args) save(simgridalignsqw1.1, file = paste0(simfolder, 'simgridalignsqw1.1.RData')) save(simgridalignsqw2.1, file = paste0(simfolder, 'simgridalignsqw2.1.RData')) save(simgridalignsqw3.1, file = paste0(simfolder, 'simgridalignsqw3.1.RData'))
# compact table square load(file = paste0(simfolder, 'simgridalignsqw1.1.RData')) load(file = paste0(simfolder, 'simgridalignsqw2.1.RData')) load(file = paste0(simfolder, 'simgridalignsqw3.1.RData')) simsum(list(Aligned = simgridalignsqw1.1, Oblique = simgridalignsqw2.1, Random = simgridalignsqw3.1))
Reduced spacing ($\sigma$ instead of $2\sigma$)
10 x 10 grid
tr <- make.grid(10,10, spacing = 25, detector = 'proximity') simrandomBVNs100.1 <- vector('list', 4) names(simrandomBVNs100.1) <- 1:4 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, theta = NULL, D = 4, CL = TRUE, detectfn = 'HHN', details = details, seed = 347, type = 'BVN', lambda = 0.4) for (i in 1:4) { args <- c(baseargs, list(sigmaX = 25/i^0.5, sigmaY = 25*i^0.5)) simrandomBVNs100.1[[i]] <- do.call(runEllipseSim, args) } save(simrandomBVNs100.1, file = paste0(simfolder, 'simrandomBVNs100.1.RData'))
Straight line
tr <- make.grid(36, 1, spacing = 25, detector = 'proximity') simgridalignlinw1.1 <- vector('list', 6) names(simgridalignlinw1.1) <- c(1:4,'1-4','1,4') simgridalignlinw4.1 <- simgridalignlinw3.1 <- simgridalignlinw2.1 <- simgridalignlinw1.1 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', theta = 0, CL = TRUE, detectfn = 'HHN', details = details) for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) # bvn perpendicular to traps args$theta <- 0; simgridalignlin1.3[[i]] <- do.call(runEllipseSim, args) # bvn oblique to traps args$theta <- pi/4; simgridalignlin2.3[[i]] <- do.call(runEllipseSim, args) # bvn parallel to traps args$theta <- pi/2; simgridalignlin3.3[[i]] <- do.call(runEllipseSim, args) # bvn random orientation args$theta <- NULL; simgridalignlin4.3[[i]] <- do.call(runEllipseSim, args) } args <- c(baseargs, list(sigmaX = rs2xy)) args$theta <- 0; simgridalignlin1.3[[5]] <- do.call(runEllipseSim, args) args$theta <- pi/4; simgridalignlin2.3[[5]] <- do.call(runEllipseSim, args) args$theta <- pi/2; simgridalignlin3.3[[5]] <- do.call(runEllipseSim, args) args$theta <- NULL; simgridalignlin4.3[[5]] <- do.call(runEllipseSim, args) args <- c(baseargs, list(sigmaX = rs2xy2)) args$theta <- 0; simgridalignlin1.3[[6]] <- do.call(runEllipseSim, args) args$theta <- pi/4; simgridalignlin2.3[[6]] <- do.call(runEllipseSim, args) args$theta <- pi/2; simgridalignlin3.3[[6]] <- do.call(runEllipseSim, args) args$theta <- NULL; simgridalignlin4.3[[6]] <- do.call(runEllipseSim, args) save(simgridalignlin1.3, file = paste0(simfolder, 'simgridalignlin1.3.RData')) save(simgridalignlin2.3, file = paste0(simfolder, 'simgridalignlin2.3.RData')) save(simgridalignlin3.3, file = paste0(simfolder, 'simgridalignlin3.3.RData')) save(simgridalignlin4.3, file = paste0(simfolder, 'simgridalignlin4.3.RData'))
Hollow square
tr <- make.grid(10, 10, spacing = 25, detector = 'proximity', hollow = TRUE) simgridalignsq1.1 <- vector('list', 4) # 'w' for wide names(simgridalignsq1.1) <- 1:4 simgridalignsq3.1 <- simgridalignsq2.1 <- simgridalignsq1.1 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = TRUE, detectfn = 'HHN', details = details) for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) # bvn aligned to traps args$theta <- 0; simgridalignsq1.1[[i]] <- do.call(runEllipseSim, args) # bvn oblique to traps args$theta <- pi/4; simgridalignsq2.1[[i]] <- do.call(runEllipseSim, args) # bvn random orientation args$theta <- NULL; simgridalignsq3.1[[i]] <- do.call(runEllipseSim, args) } args <- c(baseargs, list(sigmaX = rs2xy)) args$theta <- 0; simgridalignsq1.1[[5]] <- do.call(runEllipseSim, args) args$theta <- pi/4; simgridalignsq2.1[[5]] <- do.call(runEllipseSim, args) args$theta <- NULL; simgridalignsq3.1[[5]] <- do.call(runEllipseSim, args) args <- c(baseargs, list(sigmaX = rs2xy2)) args$theta <- 0; simgridalignsq1.1[[6]] <- do.call(runEllipseSim, args) args$theta <- pi/4; simgridalignsq2.1[[6]] <- do.call(runEllipseSim, args) args$theta <- NULL; simgridalignsq3.1[[6]] <- do.call(runEllipseSim, args) save(simgridalignsq1.2, file = paste0(simfolder, 'simgridalignsq1.2.RData')) save(simgridalignsq2.2, file = paste0(simfolder, 'simgridalignsq2.2.RData')) save(simgridalignsq3.2, file = paste0(simfolder, 'simgridalignsq3.2.RData'))
Spacing comparisons
10 x 10 grid
load(file = paste0(simfolder, 'simrandomBVN100.3.RData')) load(file = paste0(simfolder, 'simrandomBVNs100.1.RData')) par(mfrow=c(1,2)) simplot(list(BVN = simrandomBVNs100.1), legend = FALSE) mtext(side=3, text = '25 m') simplot(list(BVN = simrandomBVN100.3[1:4]), legend = FALSE) mtext(side=3, text = '50 m')
Straight line
load(file = paste0(simfolder, 'simgridalignlin1.3.RData')) load(file = paste0(simfolder, 'simgridalignlin2.3.RData')) load(file = paste0(simfolder, 'simgridalignlin3.3.RData')) load(file = paste0(simfolder, 'simgridalignlin4.3.RData')) load(file = paste0(simfolder, 'simgridalignlinw1.1.RData')) load(file = paste0(simfolder, 'simgridalignlinw2.1.RData')) load(file = paste0(simfolder, 'simgridalignlinw3.1.RData')) load(file = paste0(simfolder, 'simgridalignlinw4.1.RData')) par(mfrow=c(1,2)) simplot(list(Perpendicular = simgridalignlin1.3[1:4], Oblique = simgridalignlin2.3[1:4], Parallel = simgridalignlin3.3[1:4], Random = simgridalignlin4.3[1:4]), legend = TRUE, ylim = c(-1,3.5)) mtext(side=3, text = '25 m') simplot(list(Perpendicular = simgridalignlinw1.1[1:4], Oblique = simgridalignlinw2.1[1:4], Parallel = simgridalignlinw3.1[1:4], Random = simgridalignlinw4.1[1:4]), legend = TRUE, ylim = c(-1,3.5)) mtext(side=3, text = '50 m')
Hollow square
load(file = paste0(simfolder, 'simgridalignsq1.2.RData')) load(file = paste0(simfolder, 'simgridalignsq2.2.RData')) load(file = paste0(simfolder, 'simgridalignsq3.2.RData')) load(file = paste0(simfolder, 'simgridalignsqw1.1.RData')) load(file = paste0(simfolder, 'simgridalignsqw2.1.RData')) load(file = paste0(simfolder, 'simgridalignsqw3.1.RData')) par(mfrow =c(1,2)) simplot(list(Aligned = simgridalignsq1.2[1:4], Oblique = simgridalignsq2.2[1:4], Random = simgridalignsq3.2[1:4]), legend = TRUE, ylim = c(-1,3.5)) mtext(side=3, text = '25 m') simplot(list(Aligned = simgridalignsqw1.1[1:4], Oblique = simgridalignsqw2.1[1:4], Random = simgridalignsqw3.1[1:4]), legend = TRUE, ylim = c(-1,3.5)) mtext(side=3, text = '50 m')
Variations
Low density (0.5/ha)
tr <- make.grid(10, 10, spacing = 50, detector = 'proximity') simrandom100low.2 <- vector('list', 4) names(simrandom100low.2) <- 1:4 simrandomBVN100low.2 <- simrandom100low.2 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, D = 0.5, CL = TRUE, detectfn = 'HHN', details = details) for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) args$type <- 'uniform'; args$g0 <- 0.2 simrandom100low.2[[i]] <- do.call(runEllipseSim, args) args$type <- 'BVN'; args$lambda0 <- 0.4 simrandomBVN100low.2[[i]] <- do.call(runEllipseSim, args) } save(simrandom100low.2, file = paste0(simfolder, 'simrandom100low.2.RData')) save(simrandomBVN100low.2, file = paste0(simfolder, 'simrandomBVN100low.2.RData'))
load(file = paste0(simfolder, 'simrandom100low.2.RData')) load(file = paste0(simfolder, 'simrandomBVN100low.2.RData')) par(mfrow = c(1,1), mar = c(4,4,4,4), xpd = FALSE) simplot(list(Uniform = simrandom100low.2, BVN = simrandomBVN100low.2), trueval = 0.5, legend = TRUE) simsum(list(Uniform = simrandom100low.2, BVN = simrandomBVN100low.2), trueval = 0.5, compact = NULL)
Small array (6 x 6 grid)
tr <- make.grid(6, 6, spacing = 50, detector = 'proximity') simrandom36.2 <- vector('list', 4) names(simrandom36.2) <- 1:4 simrandomBVN36.2 <- simrandom36.2 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, D = 4, CL = TRUE, detectfn = 'HHN', details = details) for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) args$type <- 'uniform'; args$g0 <- 0.2 simrandom36.2[[i]] <- do.call(runEllipseSim, args) args$type <- 'BVN'; args$lambda0 <- 0.4 simrandomBVN36.2[[i]] <- do.call(runEllipseSim, args) } save(simrandom36.2, file = paste0(simfolder, 'simrandom36.2.RData')) save(simrandomBVN36.2, file = paste0(simfolder, 'simrandomBVN36.2.RData'))
load(file = paste0(simfolder, 'simrandom36.2.RData')) load(file = paste0(simfolder, 'simrandomBVN36.2.RData')) par(mfrow = c(1,1), mar = c(4,4,4,4), xpd = FALSE) simplot(list(Uniform = simrandom36.2[1:4], BVN = simrandomBVN36.2[1:4]), legend = TRUE) simsum(list(Uniform = simrandom36.2, BVN = simrandomBVN36.2))
Common random orientation
tr <- make.grid(10, 10, spacing = 50, detector = 'proximity') simrandom100C.2 <- vector('list', 4) names(simrandom100C.2) <- 1:4 simrandomBVN100C.2 <- simrandom100C.2 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, theta = -1, D = 4, CL = TRUE, detectfn = 'HHN', details = details) for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) args$type <- 'uniform'; args$g0 <- 0.2 simrandom100C.2[[i]] <- do.call(runEllipseSim, args) args$type <- 'BVN'; args$lambda0 <- 0.4 simrandomBVN100C.2[[i]] <- do.call(runEllipseSim, args) } save(simrandom100C.2, file = paste0(simfolder, 'simrandom100C.2.RData')) save(simrandomBVN100C.2, file = paste0(simfolder, 'simrandomBVN100C.2.RData'))
load(file = paste0(simfolder, 'simrandom100C.2.RData')) load(file = paste0(simfolder, 'simrandomBVN100C.2.RData')) simplot(list(Uniform = simrandom100C.2, BVN = simrandomBVN100C.2), legend = TRUE)
Grid with varying orientation
Spacings narrow and wide (w); orientation parallel (p) and oblique(o). Only consider aspect ratio 4 for now.
tr <- make.grid(10,10, spacing = 25, detector = 'proximity') simrandomBVN100p.1 <- vector('list', 4) names(simrandomBVN100p.1) <- 1:4 simrandomBVN100o.1 <- simrandomBVN100p.1 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, type = 'BVN', lambda0 = 0.4, D = 4, CL = TRUE, detectfn = 'HHN', details = details) for (i in 4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY, theta = 0)) simrandomBVN100p.1[[i]] <- do.call(runEllipseSim, args) args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY, theta = pi/4)) simrandomBVN100o.1[[i]] <- do.call(runEllipseSim, args) } tr <- make.grid(10,10, spacing = 50, detector = 'proximity') simrandomBVN100pw.1 <- vector('list', 4) names(simrandomBVN100pw.1) <- 1:4 simrandomBVN100ow.1 <- simrandomBVN100pw.1 for (i in 4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY, theta = 0)) simrandomBVN100pw.1[[i]] <- do.call(runEllipseSim, args) args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY, theta = pi/4)) simrandomBVN100ow.1[[i]] <- do.call(runEllipseSim, args) } save(simrandomBVN100p.1, file = paste0(simfolder, 'simrandomBVN100p.1.RData')) save(simrandomBVN100o.1, file = paste0(simfolder, 'simrandomBVN100o.1.RData')) save(simrandomBVN100pw.1, file = paste0(simfolder, 'simrandomBVN100pw.1.RData')) save(simrandomBVN100ow.1, file = paste0(simfolder, 'simrandomBVN100ow.1.RData'))
Anisotropic model
Using function anisotropic.fit. Extract the fitted coefficients corresponding to psiA and psiR with the coef method for secr objects, and rely on direct estimation of density (CL = FALSE).
Hollow square
tr <- make.grid(10, 10, spacing = 25, detector = 'proximity', hollow = TRUE) simaniso1 <- vector('list', 4); names(simaniso1) <- 1:4 simaniso4 <- simaniso3 <- simaniso2 <- simaniso1 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN', details = list(distribution = 'binomial'), extractfn = coef, secrfn = 'anisotropic.fit') for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) args$theta <- 0 # aligned to detectors simaniso1[[i]] <- do.call(runEllipseSim, args) args$theta <- pi/4 # oblique to detectors simaniso2[[i]] <- do.call(runEllipseSim, args) args$theta <- NULL # random orientation simaniso3[[i]] <- do.call(runEllipseSim, args) args$theta <- -1 # common random orientation simaniso4[[i]] <- do.call(runEllipseSim, args) } save(simaniso1, file = paste0(simfolder, 'simaniso1.RData')) save(simaniso2, file = paste0(simfolder, 'simaniso2.RData')) save(simaniso3, file = paste0(simfolder, 'simaniso3.RData')) save(simaniso4, file = paste0(simfolder, 'simaniso4.RData'))
load(file = paste0(simfolder, 'simgridalignsq1.2.RData')) load(file = paste0(simfolder, 'simgridalignsq2.2.RData')) load(file = paste0(simfolder, 'simgridalignsq3.2.RData')) load(file = paste0(simfolder, 'simaniso1.RData')) load(file = paste0(simfolder, 'simaniso2.RData')) load(file = paste0(simfolder, 'simaniso3.RData')) load(file = paste0(simfolder, 'simaniso4.RData')) par(mfrow = c(1,2)) simplot(list(Aligned = simgridalignsq1.2[1:4], Oblique = simgridalignsq2.2[1:4], Random = simgridalignsq3.2[1:4]), ylim = c(-0.4,0.45), legend = TRUE, pchi = c(23, 22, 16)) text(-0.4, 0.5, 'a.', cex = 1.45, xpd = TRUE) text(3, -0.35, 'isotropic') simplot(list(Aligned = simaniso1, Oblique = simaniso2, Random = simaniso3), component = "pred", ylim = c(-0.4,0.45), legend = TRUE, pchi = c(23, 22, 16)) text(-0.4, 0.5, 'b.', cex = 1.45, xpd = TRUE) text(3, -0.35, 'anisotropic')
Tabular summary for anisotropic model.
simsum(list(Aligned = simaniso1, Oblique = simaniso2, Random = simaniso3), component = "pred")
Tabular comparison of isotropic and anisotropic models for hollow grid.
tab1 <- simsum(list(Aligned = simgridalignsq1.2[1:4], Oblique = simgridalignsq2.2[1:4], Random = simgridalignsq3.2[1:4])) # find D-hat in the 'fit' component tab2 <- simsum(list(Aligned = simaniso1, Oblique = simaniso2, Random = simaniso3), component = "pred") # find D-hat in the 'pred' component fn <- function(t1,t2) cbind(t1[,-1],t2[,-1]) mapply(fn, tab1, tab2, SIMPLIFY = FALSE)
Straight line
tr <- make.grid(36, 1, spacing = 25, detector = 'proximity') simanisolin1 <- vector('list', 4); names(simanisolin1) <- 1:4 simanisolin5 <- simanisolin4 <- simanisolin3 <- simanisolin2 <- simanisolin1 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN', details = list(distribution = 'binomial'), extractfn = coef, secrfn = 'anisotropic.fit') for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) args$theta <- 0 # perpendicular to detectors simanisolin1[[i]] <- do.call(runEllipseSim, args) args$theta <- pi/4 # oblique to detectors simanisolin2[[i]] <- do.call(runEllipseSim, args) args$theta <- NULL # random orientation simanisolin3[[i]] <- do.call(runEllipseSim, args) args$theta <- -1 # common random orientation simanisolin4[[i]] <- do.call(runEllipseSim, args) args$theta <- pi/2 # parallel simanisolin5[[i]] <- do.call(runEllipseSim, args) # partial run 2018-08-23 } save(simanisolin1, file = paste0(simfolder, 'simanisolin1.RData')) save(simanisolin2, file = paste0(simfolder, 'simanisolin2.RData')) save(simanisolin3, file = paste0(simfolder, 'simanisolin3.RData')) save(simanisolin4, file = paste0(simfolder, 'simanisolin4.RData')) save(simanisolin5, file = paste0(simfolder, 'simanisolin5.RData'))
load(file = paste0(simfolder, 'simgridalignlin1.3.RData')) load(file = paste0(simfolder, 'simgridalignlin2.3.RData')) load(file = paste0(simfolder, 'simgridalignlin4.3.RData')) load(file = paste0(simfolder, 'simanisolin1.RData')) load(file = paste0(simfolder, 'simanisolin2.RData')) load(file = paste0(simfolder, 'simanisolin3.RData')) par(mfrow = c(1,2)) simplot(list(Aligned = simgridalignlin1.3[1:4], Oblique = simgridalignlin2.3[1:4], Random = simgridalignlin4.3[1:4]), ylim = c(-1,3.5), legend = TRUE, pchi = c(23, 22, 16)) text(-0.1, 3.5, 'a.', cex = 1.45, xpd = TRUE) text(1, -0.8, 'isotropic', adj=0) simplot(list(Perpendicular = simanisolin1, Oblique = simanisolin2, Random = simanisolin3), component = "pred", ylim = c(-1,3.5), legend = TRUE, pchi = c(23, 22, 16)) text(-0.1, 3.5, 'b.', cex = 1.45, xpd = TRUE) text(1, -0.8, 'anisotropic', adj=0)
Conclude that anisotropic transformation largely ineffective with linear array.
10 x 10 grid
tr <- make.grid(10, 10, spacing = 25, detector = 'proximity') simanisogrid1 <- vector('list', 4); names(simanisogrid1) <- 1:4 simanisogrid4 <- simanisogrid3 <- simanisogrid2 <- simanisogrid1 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN', details = list(distribution = 'binomial'), extractfn = coef, secrfn = 'anisotropic.fit') for (i in 1:4) { sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 args <- c(baseargs, list(sigmaX = sigmaX, sigmaY = sigmaY)) args$theta <- 0 # aligned to detectors simanisogrid1[[i]] <- do.call(runEllipseSim, args) args$theta <- pi/4 # oblique to detectors simanisogrid2[[i]] <- do.call(runEllipseSim, args) args$theta <- NULL # random orientation simanisogrid3[[i]] <- do.call(runEllipseSim, args) args$theta <- -1 # common random orientation simanisogrid4[[i]] <- do.call(runEllipseSim, args) } save(simanisogrid1, file = paste0(simfolder, 'simanisogrid1.RData')) save(simanisogrid2, file = paste0(simfolder, 'simanisogrid2.RData')) save(simanisogrid3, file = paste0(simfolder, 'simanisogrid3.RData')) save(simanisogrid4, file = paste0(simfolder, 'simanisogrid4.RData'))
load(file = paste0(simfolder, 'simanisogrid3.RData')) load(file = paste0(simfolder, 'simrandomBVNs100.1.RData')) simsum(list(Simple = simrandomBVNs100.1)) simsum(list(Anisotropic = simanisogrid3), component = "pred")
Performs well with full grid data (narrow spacing).
Now run all geometries with wide spacing (w) for aspect ratio 4.0.
tr <- make.grid(10, 10, spacing = 50, detector = 'proximity') simanisogridw1 <- vector('list', 4); names(simanisogridw1) <- 1:4 simanisogridw3 <- simanisogridw2 <- simanisogridw1 i <- 4 sigmaX <- 25/i^0.5; sigmaY <- 25*i^0.5 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN', details = list(distribution = 'binomial'), extractfn = coef, secrfn = 'anisotropic.fit', sigmaX = sigmaX, sigmaY = sigmaY) args <- c(baseargs, list(theta = 0)) # aligned to detectors simanisogridw1[[i]] <- do.call(runEllipseSim, args) args$theta <- pi/4 # oblique to detectors simanisogridw2[[i]] <- do.call(runEllipseSim, args) args$theta <- NULL # random orientation simanisogridw3[[i]] <- do.call(runEllipseSim, args) save(simanisogridw1, file = paste0(simfolder, 'simanisogridw1.RData')) save(simanisogridw2, file = paste0(simfolder, 'simanisogridw2.RData')) save(simanisogridw3, file = paste0(simfolder, 'simanisogridw3.RData')) tr <- make.grid(36, 1, spacing = 50, detector = 'proximity') baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN', details = list(distribution = 'binomial'), extractfn = coef, secrfn = 'anisotropic.fit', sigmaX = sigmaX, sigmaY = sigmaY) simanisolinw1 <- vector('list', 4); names(simanisolinw1) <- 1:4 simanisolinw4 <- simanisolinw3 <- simanisolinw2 <- simanisolinw1 args <- c(baseargs, list(theta = 0)) # perpendicular to detectors simanisolinw1[[i]] <- do.call(runEllipseSim, args) args$theta <- pi/4 # oblique to detectors simanisolinw2[[i]] <- do.call(runEllipseSim, args) args$theta <- pi/2 # parallel to detectors simanisolinw3[[i]] <- do.call(runEllipseSim, args) args$theta <- NULL # random orientation simanisolinw4[[i]] <- do.call(runEllipseSim, args) save(simanisolinw1, file = paste0(simfolder, 'simanisolinw1.RData')) save(simanisolinw2, file = paste0(simfolder, 'simanisolinw2.RData')) save(simanisolinw3, file = paste0(simfolder, 'simanisolinw3.RData')) save(simanisolinw4, file = paste0(simfolder, 'simanisolinw4.RData')) tr <- make.grid(10, 10, spacing = 50, detector = 'proximity', hollow = TRUE) simanisow1 <- vector('list', 4); names(simanisow1) <- 1:4 simanisow3 <- simanisow2 <- simanisow1 baseargs <- list(nrepl = nrepl, buffer = 200, ncores = ncores, traps = tr, lambda0 = 0.4, D = 4, type = 'BVN', CL = FALSE, detectfn = 'HHN', details = list(distribution = 'binomial'), extractfn = coef, secrfn = 'anisotropic.fit', sigmaX = sigmaX, sigmaY = sigmaY) args <- c(baseargs, list(theta = 0)) # aligned to detectors simanisow1[[i]] <- do.call(runEllipseSim, args) args$theta <- pi/4 # oblique to detectors simanisow2[[i]] <- do.call(runEllipseSim, args) args$theta <- NULL # random orientation simanisow3[[i]] <- do.call(runEllipseSim, args) save(simanisow1, file = paste0(simfolder, 'simanisow1.RData')) save(simanisow2, file = paste0(simfolder, 'simanisow2.RData')) save(simanisow3, file = paste0(simfolder, 'simanisow3.RData'))
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.