www.phidot.org

Hi all,

I have a situation where I've calculated SIFs from model averaged occupancy outputs from a few different sets of two-species, single season models. I am having trouble figuring out how to calculate the SEs for these SIFs. In looking into sources, I found a paper by Robinson, et. al, from 2014 where they cite the delta method as how they calculating their SIF's SE. After reading through the appendix in the Mark Manual covering the method I am still having trouble figuring out how to apply it to SIFs.

Does anyone have any experience with this sort of calculation or any pointers on where I could find an example?

If you have the model SIF's and variances, from each model (should be in the model$derived object), you can compute the model-averaged SIF and std.error as:

Code: Select all

ma.SIF = sum(wgt * sif)

ma.SIF.var = sum(wgt * sif.var) + sum(wgt*(sif - ma.sif)^2)
ma.SIF.se = sqrt(ma.SIF.var)

The variance is the sum of the within model variances, sif.var, and the between model variances, (sif-ma.sif)^2.

wgt is a vector of model weights
sif is a vector of SIF estimates for each model
sif.var is a vector of var(SIF) for each model

Hi Jim,

Thanks for the lightning fast response. That all makes sense, but maybe I have something set up incorrectly because I can't seem to find the SIF output. When I check the $derived portion of my model output, there are only nu and rho values listed. None of the values listed make sense given the SIF formula and the outputs for PsiA, PsiBA, and PsiBa. From playing around in Presence and reading elsewhere on here, I think SIF is generally given as phi, right?

For all of my candidate models I am using the occMod function with the type set to "so.2sp.1" and the parameterization set to "PsiBA". Should I have them set up differently?

Thanks for any pointers you can give me,
Chris

Hi Chris
There are multiple definitions for the SIF, which are you using?

If it's the ratio of probabilities (ie PsiAB/(PsiA*PsiB) then I recommend against using that one, even though it's one that we first suggested. The nu and rho SIF's are based on the ratios of odds (instead of ratios of probabilities) and have better properties. I generally recommend people use these instead.

Cheers
Darryl

Hi Darryl,

Thanks for the help.

I was calculating my own SIFs from the PsiA, PsiBA, and PsiBa outputs using the formula recommended in Robinson et al., 2014 from Richmond et al., 2010:

SIF = (PsiA*PsiBA) / PsiA * (PsiA*PsiBA + (1 - PsiA) * PsiBa)

If this is no longer the recommended method, then I'll definitely go the route of using and reporting nu, thanks.

One question on nu, though: I have 108 sites in my sample and for each of my PsiA, PsiBA, and PsiBa outputs at these sites I am getting a different value based on the gradients of values of the covariates involved in the model. When I go to the derived output of each model, the nu values are all the same. Any idea why that might be the case?

Thanks again,
Chris

You must have defined the model where the interaction (on the odds scale) between species is constant (ie doesn't vary with any covariates). Hard to some more without knowing the details of what you've done already.

Darryl,

Turns out I had some of models defined incorrectly, just as you said. Thanks for clearing that up!

-Chris

Hi Chris,

Are you using the latest version of RPresence? I updated the output a while back for the two -species model so it outputs psiA, psiB, psiAB and phi as derived parameters.

Hi Jim,

An update helped quite a bit. Thanks to you both for the input and pointers!

- Chris

Hello, I also need some helps in estimating the CI of SIF. I'm using the "psiBA" parameterization, and calculating SIF through this formula (Richmond et al. 2010) as well:

SIF = (PsiA*PsiBA) / PsiA * (PsiA*PsiBA + (1 - PsiA) * PsiBa)

jhines wrote:If you have the model SIF's and variances, from each model (should be in the model$derived object), you can compute the model-averaged SIF and std.error as:

Code: Select all

ma.SIF = sum(wgt * sif)

ma.SIF.var = sum(wgt * sif.var) + sum(wgt*(sif - ma.sif)^2)
ma.SIF.se = sqrt(ma.SIF.var)

The variance is the sum of the within model variances, sif.var, and the between model variances, (sif-ma.sif)^2.

wgt is a vector of model weights
sif is a vector of SIF estimates for each model
sif.var is a vector of var(SIF) for each model

The formula you provided to calculate model-averaged SIF really helps. However, I am not sure if I understand the sif.var correctly.

Code: Select all

fit <- occMod(model=list(psi ~ SP*elevation+INT*roaddensity,
                            p ~ SP+INT_o+INT_d+SP:INT_o),
                 data = dog_mongoose.pao, type="so.2sp.1", parameterization="PsiBA")

psiBa <- as.vector(fit$real$psi$psiBa[,1])
psiBA <- as.vector(fit$real$psi$psiBA[,1])
psiA <- as.vector(fit$real$psi$psiA[,1])

#a vector containing SIF for each camera
SIF = (psiA*psiBA) / psiA * (psiA*psiBA + (1 - psiA) * psiBa) 
sif.var = var(SIF)
sif.se = sqrt(sif.var)
sif = mean(SIF)  # sif for this model
mean(SIF+1.96*sif.se) # upper CI
mean(SIF-1.96*sif.se) # lower CI

The fit$derived output only contains psiA, psiB, psiAB, nu, and phi. I can not find psiBa in the output, so I store the psiBa for each camera as "as.vector(fit$real$psi$psiBa[,1])" instead.
I first calculate the SIFs for every camera (the vector "SIF"), and I regard the variance within the vector SIF as sif.var of this model. Is that correct?

Thank you for reading!