Alternate Parts and Density Covariates

[Clayton Lamb]

Hello,

I am interested in linking sigma with density. I have a large dataset (2,500 animals) covering a large density gradient and many sessions (49). The sigma and a0 or c methods detailed in https://www.otago.ac.nz/density/pdfs/secr-parameterisations.pdf seem like a good fit for  my work.

I am curious why this only works when D~session, but not when D~mask covariate. If, for example, D~mask covariate produces a better model, shouldn’t this produce better density, and sigma estimates for each session? Clearly there is some structural reason why this doesn’t work at the moment, but Im not sure why.

EX:

library(secr)

##make mask

clippedmask <- make.mask(traps(ovenCH), type = 'trapbuffer', buffer = 300, spacing=20)

##add fake mask covariates

for (i in 1:length(clippedmask)){

  covariates(clippedmask[[i]]) <- data.frame(rand = runif(1850, -2.0, 2))

}

##model from https://www.otago.ac.nz/density/pdfs/secr-parameterisations.pdf

oven0509bs <- secr.fit(ovenCH, mask=clippedmask, model = list(D ~ session, sigmak ~ 1,

                                            a0 ~ bk), buffer = 300, detectfn = 'HHN', trace = FALSE)

##with density covariate

oven0509dc <- secr.fit(ovenCH, mask=clippedmask, model = list(D ~ rand, sigmak ~ 1,

                                                              a0 ~ bk), buffer = 300, detectfn = 'HHN', trace = FALSE)

###Error in secr.fit(ovenCH, mask = clippedmask, model = list(D ~ rand, sigmak ~  : 

##only session and group models allowed for density when details$param = 5

Thanks!

--
Clayton Lamb, PhD Candidate
South Rockies Grizzly Bear Project
S.Boutin and S.Nielsen Labs
Biological Sciences, University of Alberta 
e: ctl...@ualberta.ca
p: 778-215-0334
t: @ClaytonTLamb
w: www.lamb-eco-research.ca

[Murray Efford]

I'm not remembering all the details, but -

What you are suggesting implicitly makes sigma a function of the mask covariates, when sigma is generally fixed across the mask. We used a workaround for spatially determined density-dependent sigma in a possum example in the Ecography paper - using x and y as predictors (not sure about other covariates, but I think they work). See Appendix 3 in the online appendix.

Murray

--
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.
For more options, visit https://groups.google.com/d/optout.

[Clayton Lamb]

Thanks, Murray. This approach in Apx 3 seems to be what I am looking for, thanks!

I used both approaches suggested in A3 ( opt1: along least cost path, and opt2: dens-dep sigma at endpoints), but both ended up throwing errors for my two continuous mask covariates (they run fine normally):

OPT2

 rd.ndvi <- secr.fit(grizzCH, model = list(D~dndvi + roaddens, noneuc ~ 1), mask = GBmask2000b, ncores=4, detectfn = 'HHN', details =list(userdist = userdist2), link =list(noneuc = 'log'), fixed=list(sigma=1))

Checking data 

Error in secr.fit(grizzCH, model = list(D ~ dndvi + roaddens, noneuc ~  : 

  parameters in model not consistent with detectfn etc. : noneuc

OPT1

rd.ndvi <- secr.fit(grizzCH, model = list(D~dndvi + roaddens, noneuc ~ 1), mask = GBmask2000b, ncores=4, detectfn = 'HHN', details =list(userdist = userdist1), link =list(noneuc = 'log'), fixed=list(sigma=1))

Checking data 

Preparing detection design matrices 

Preparing density design matrix 

Preparing non-Euclidean parameter design matrix 

Finding initial parameter values... 

Initial values  D = 0.00027, g0 = 0.28788, sigma = 2943.91937 

Maximizing likelihood... 

Eval     Loglik          D    D.dndvi D.roaddens    lambda0     noneuc 

Error in checkForRemoteErrors(val) : 

  4 nodes produced errors; first error: could not find function “covariates"

Clayton

[Murray Efford]

Suggest you forget about using ncores>1 (it looks like that crashes OPT1 - may be fixable, but I don't have time).

In OPT2 I wonder if your userdist2 correctly identifies the parameters (see secr-noneuclidean.pdf 'Dynamic userdist').

To unsubscribe from this group and stop receiving emails from it, send an email to secrgroup+unsubscribe@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

[Clayton Lamb]

Awesome, thanks Murray. Removing the multi-cores fixed the problem for OPT1. I’ll paste a working, reproducible script of this below, which builds on the OvenCH example I used to start this thread.

I wasn’t able to get Opt2 working. 

library(secr)

##make mask

clippedmask <- make.mask(traps(ovenCH), type = 'trapbuffer', buffer = 300, spacing=20)

##add fake mask covariates

for (i in 1:length(clippedmask)){

  covariates(clippedmask[[i]]) <- data.frame(rand = runif(1850, -2.0, 2))

}

##model from https://www.otago.ac.nz/density/pdfs/secr-parameterisations.pdf

oven0509bs <- secr.fit(ovenCH, mask=clippedmask, model = list(D ~ session, sigmak ~ 1,

                                            a0 ~ bk), buffer = 300, detectfn = 'HHN')

###WORKING EXAMPLE TO LINK SIGMA WITH DENSITY PREDICTED FROM MASK COVARIATES ALONE

##Define function for non-euclidian distance

##from http://www.ecography.org/sites/ecography.org/files/appendix/ecog-01511.pdf

## Option 1: integrating along least-cost path

userdist1 <- function (xy1, xy2, mask) {

  if (missing(xy1)) return(c("D", "noneuc"))

  D <- covariates(mask)$D

  k <- covariates(mask)$noneuc ## noneuc = k

  sigma <- k / sqrt(D)

  if (!require(gdistance))

    stop ('install package gdistance to use this function')

  ## make raster from mask, assigning 1/sigma to all habitat cells

  Draster <- raster(mask, values = 1/sigma)

  ## Assume animals can traverse gaps: bridge gaps using mean

  Draster[is.na(Draster[])] <- mean(Draster[], na.rm = TRUE)

  ## form gdistance ‘TransitionLayer’ to define cell-cell connections

  trfn <- function(x) 1 / mean(x)

  tr <- transition(Draster, transitionFunction = trfn, directions = 16)

  tr <- geoCorrection(tr, type = "c", multpl = FALSE)

  ## costDistance finds length of least-cost path

  distmat <- costDistance(tr, as.matrix(xy1), as.matrix(xy2))

  distmat[!is.finite(distmat)] <- 1e10 ## avoid infinite values

  distmat

}

o1 <- secr.fit(ovenCH, mask=clippedmask,model = list(D~rand, noneuc ~ 1), detectfn = 'HHN', details =list(userdist = userdist1), link =

                      list(noneuc = 'log'), fixed=list(sigma=1))

## plotting

addSigma <- function (Dsurface, Dcovariate = 'D.0', k) {

  D <- covariates(Dsurface)[,Dcovariate]

  covariates(Dsurface)[,'sigma'] <- k / sqrt(D)

  Dsurface

}

DS <- predictDsurface(o1)

DS <- addSigma(DS[[1]], k = predict(o1)[[1]]$estimate[3] )

plot(covariates(DS)$D.0, covariates(DS)$sigma)

###NOT WORKING YET

## Option 2: averaging density-dependent sigma at endpoints of path

userdist2 <- function (xy1, xy2, mask) {

  ## use addCovariates here to associate D with each point

  xy1 <- addCovariates(xy1, mask)

  xy2 <- addCovariates(xy2, mask)

  D1 <- covariates(xy1)$D^0.5

  D2 <- covariates(xy2)$D^0.5

  k <- covariates(mask)$noneuc[1] ## assume k constant

  euc <- edist(xy1, xy2)

  sig <- outer (D1, D2, FUN=function(d1,d2) k * (1/d1 + 1/d2)/2)

  euc / sig

}

o2 <- secr.fit(ovenCH, mask=clippedmask,model = list(D~rand, noneuc ~ 1), detectfn = 'HHN', details =list(userdist = userdist2), link =

                 list(noneuc = 'log'), fixed=list(sigma=1))

To unsubscribe from this group and stop receiving emails from it, send an email to secrgroup+...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

[Clayton Lamb]

Just following up on this, I have applied this approach (OPT 1) to my data. The models threw the “at least one variance calculation failed error”, but I resolved by re-running these models using method=“none”.  So all seemed good. When I print the model, I get an odd error that I don’t get on the multi-session worked example provided earlier in this thread. I can’t seem to track down the issue. I am able to predict a density surface with this model, so not sure why the summary is failing to print. Just want to make sure this doesn’t signify an underlying issue in my model

Fitted (real) parameters evaluated at base levels of covariates 

Error in lpred[, 1] <- Yp : 

  number of items to replace is not a multiple of replacement length

In addition: Warning message:

In secr.lpredictor(formula = models[[x]], newdata = newdata, indx = parindices[[x]],  :

  missing values in predictors?

Cheers,

Clayton

<Screen Shot 2019-01-27 at 11.16.10 AM.png>