Trap covariates in simulations

[stefano palmero]

Hi all,

I would like to include trap types (on/off road) in my SECRDESIGN scenarios but I am not sure how to do it.

First of all, trap types should be specified in the trap file.

When I know type-specific g0s, should I create two groups? Something like this:

## Scenarios
densities <- c(5e-05,1e-04,1.5e-04,2e-04)
g0 <- 0.01
sigmas <- c(3000,4000,5000)

scen <- make.scenarios(noccasions = c(60,80,100,120), nrepeats = 1,
                       D = densities * 0.75,
                       sigma = sigmas,
                       g0 = g0,
                       groups = c('on','off'))
off <- scen$group == 'off'
scen$g0[off] <- 0.02

and then of course g0~traptype etc...

But I am not sure about communication between trap file, scenarios and the model in this case

Thank you very much

Stefano

[Murray Efford]

Hi Stefano

Groups could work (examples at https://github.com/MurrayEfford/secr-simulations under STR). But there is a more direct way using explicit input of detector-specific detection parameters in sim.capthist. For a different reason (aiming to reduce number of occasions)  I started out with HHN instead of HN (hence lambda0 instead of g0), but this should also work with HN/g0. It's not pretty, but it seems to get there.

Murray

library(secrdesign)

# some detectors with a binary covariate 1,2
grid<- make.grid(8,8, spacing = 5000, detector = 'proximity')
onoff <- sample.int(2, nrow(grid), replace = TRUE)
covariates(grid) <- data.frame(onoff = onoff)
plot(secr::as.mask(grid), cov='onoff', dots = FALSE)

lambdas <- c(0.01,0.02)  # for covariate levels 1,2
sigmas  <- c(3000,4000,5000)

# manually construct some occasion x detector matrix inputs for sim.capthist

# (ultimately replaces lambda0, sigma from scen)
onedetpar <- function (noccasions, sigma) {
# filling by row re-uses detector-specific covariates on each occasion
mat <- matrix(lambdas[onoff], byrow = TRUE, nrow = noccasions, ncol = nrow(grid))
    list(detectpar = list(lambda0 = mat, sigma = sigma))
}

# all combinations of duration and sigma

dets <- expand.grid(noccasions = c(60,80,100,120), sigma = sigmas)
detarglist <- mapply(onedetpar, dets$noccasions, dets$sigma, SIMPLIFY = FALSE)

## Scenarios
densities <- c(5e-05,1e-04,1.5e-04,2e-04)

scen <- make.scenarios(
noccasions = 60,       # placeholder - replaced later

D          = densities * 0.75,

detectfn   = 'HHN',
lambda0    = 0.01,     # placeholder - replaced later
sigma      = 3000,     # placeholder - replaced later
detindex   = 1:12)

# force noccasions to follow matrices in detarg
scen$noccasions <- dets$noccasions[scen$detindex]

sims <- run.scenarios(
nrepl     = 20,
scenarios = scen,
trapset   = grid,
det.args  = detarglist,
fit       = TRUE,
fit.args  = list(detectfn = 'HHN', model = lambda0~onoff))

estimateSummary(sims)

--
You received this message because you are subscribed to the Google Groups "secr" group.
To unsubscribe from this group and stop receiving emails from it, send an email to secrgroup+...@googlegroups.com.
To view this discussion, visit https://groups.google.com/d/msgid/secrgroup/8310ff0e-f1c7-42a6-ad83-d649e57f4a48n%40googlegroups.com.

[stefano palmero]

Dear Murray,

Thank you very much for the help, as always.

I'll try this code out

Best

Stefano Palmero

Faculty of Environment and Natural Resources

University of Freiburg

Tennenbacher Straße 4

79106 Freiburg

Germany

[stefano palmero]

Hi Murray,

I finally tried out this code and it works.

Only I cannot see the results for both levels of my covariates for g0.

For example, when running estimateSummary(sims), I can see both levels in the det.args but not in the estimates:

$det.args
$det.args[[1]]
$det.args[[1]]$detectpar
$det.args[[1]]$detectpar$g0
     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18]
[1,]  0.2  0.2 0.07 0.07 0.07 0.07  0.2 0.07 0.07  0.07  0.07   0.2  0.07  0.07  0.07   0.2  0.07  0.07
[2,]  0.2  0.2 0.07 0.07 0.07 0.07  0.2 0.07 0.07  0.07  0.07   0.2  0.07  0.07  0.07   0.2  0.07  0.07
[3,]  0.2  0.2 0.07 0.07 0.07 0.07  0.2 0.07 0.07  0.07  0.07   0.2  0.07  0.07  0.07   0.2  0.07  0.07
[4,]  0.2  0.2 0.07 0.07 0.07 0.07  0.2 0.07 0.07  0.07  0.07   0.2  0.07  0.07  0.07   0.2  0.07  0.07
[5,]  0.2  0.2 0.07 0.07 0.07 0.07  0.2 0.07 0.07  0.07  0.07   0.2  0.07  0.07  0.07   0.2  0.07  0.07
[6,]  0.2  0.2 0.07 0.07 0.07 0.07  0.2 0.07 0.07  0.07  0.07   0.2  0.07  0.07  0.07   0.2  0.07  0.07
     [,19] [,20] [,21] [,22] [,23] [,24] [,25] [,26] [,27] [,28] [,29] [,30] [,31] [,32] [,33]
[1,]  0.07  0.07  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07
[2,]  0.07  0.07  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07
[3,]  0.07  0.07  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07
[4,]  0.07  0.07  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07
[5,]  0.07  0.07  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07
[6,]  0.07  0.07  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07  0.07   0.2  0.07  0.07  0.07

$output

$output$`1`
$output$`1`[[1]]
       link     estimate  SE.estimate          lcl          ucl
D       log 5.272753e-02  0.001991611 4.896633e-02 5.677764e-02
g0    logit 7.126595e-02  0.003156678 6.532154e-02 7.770635e-02
sigma   log 1.938706e+03 41.810339884 1.858476e+03 2.022400e+03

I would use predict if it was a secr.fit output, but how can I make sure that the covariate was actually used?

Thank you very much

Dr. Stefano Palmero

Department of National Park Monitoring and Animal Management 

Freyunger Str. 2

94481 Grafenau

Germany

[Murray Efford]

Hi Stefano

You can specify a different extractfn, such as coef. Then the individual coefficients (beta parameters) are captured and something like  estimateSummary(sims2, 'lambda0.onoff') is informative. Try it out with small nrepl to see how it goes.

Murray