distsamp population estimates not changing with different covariates

[Christina Hunt]

Hi everyone,

I'm relatively new to distance sampling and am still getting to grips with basic models in distsamp. I've come across some unusual results and I'm not sure whether they're caused by my data or something in my code. I've tried searching in this group and elsewhere but not found anything useful to help me.

In a simple model with or without detection covariates, when I vary the lambda covariates the total population estimate calculated using sum(predict(model, type="state", newdata=NULL)[,1]) is the same, or varies by less than 1 (e.g 120.9, 121.0, 121.6). I checked the estimates for each of my sites and the values do vary between models, but weirdly they still sum to the same value. Could there be something I'm doing that is constraining the total population estimate, such that even though site estimates differ they still sum to the same value?

I have 201 sites, each sampled once. There are 50 individuals detected. This is the kind of model I'm running: 

distsamp(~1~April_temp_degrees_celsius,umf, keyfun="halfnorm",output="abund")

Any suggestions or advice would be greatly appreciated! :)

[Ken Kellner]

Do the lambda covariates appear to actually be informative in predicting abundance? Can you show the summary output of a couple of the models?

Ken

[Christina Hunt]

Hi Ken,

Yes they do appear to be informative. I have copied the code below for a couple of univariate models and for a model with both covariates in. You can see that although the site estimates change, the sum of the site estimates stays constant between models. I have also included the plotEffects graphs for the univariate models. Hopefully this may be enough to diagnose the problem. 

> m1<-distsamp(~1~April_temp_degrees_celsius,umf, keyfun="halfnorm",output="abund") > sum(predict(m1, type="state", newdata=NULL)[,1]) [1] 92.00112 > predict(m1, type="state", newdata=NULL)[1:10,1] [1] 2.4653172 2.8637173 0.9111429 0.9086730 0.9604130 1.1930856 0.6446510 0.3903778 0.5498071 0.4265304 > summary(m1) Call: distsamp(formula = ~1 ~ April_temp_degrees_celsius, data = umf, keyfun = "halfnorm", output = "abund") Abundance (log-scale): Estimate SE z P(>|z|) (Intercept) -0.959 0.201 -4.78 1.74e-06 April_temp_degrees_celsius -0.591 0.119 -4.99 6.13e-07 Detection (log-scale): Estimate SE z P(>|z|) 3.28 0.138 23.8 2.88e-125 AIC: 394.1063 Number of sites: 201 Survey design: line-transect Detection function: halfnorm UnitsIn: m UnitsOut: ha

> m2<-distsamp(~1~Percentage_forest,umf, keyfun="halfnorm",output="abund") > sum(predict(m2, type="state", newdata=NULL)[,1]) [1] 92.00114 > predict(m2, type="state", newdata=NULL)[1:10,1] [1] 0.2316300 0.2391444 1.2599119 0.2337286 0.2316300 0.2316300 0.2316300 0.2316300 0.2316300 0.2316300 > summary(m2) Call: distsamp(formula = ~1 ~ Percentage_forest, data = umf, keyfun = "halfnorm", output = "abund") Abundance (log-scale): Estimate SE z P(>|z|) (Intercept) -1.037 0.208 -4.97 6.60e-07 Percentage_forest 0.626 0.113 5.53 3.18e-08 Detection (log-scale): Estimate SE z P(>|z|) 3.28 0.138 23.8 2.88e-125 AIC: 387.5083 Number of sites: 201 Survey design: line-transect Detection function: halfnorm UnitsIn: m UnitsOut: ha

> m3<-distsamp(~1~April_temp_degrees_celsius+Percentage_forest,umf, keyfun="halfnorm",output="abund") > sum(predict(m3, type="state", newdata=NULL)[,1]) [1] 92.00113 > predict(m3, type="state", newdata=NULL)[1:10,1] [1] 1.3786630 1.6588944 2.3119748 0.4923582 0.5172517 0.6481506 0.3417103 0.2028260 0.2895910 0.2223944 > summary(m3) Call: distsamp(formula = ~1 ~ April_temp_degrees_celsius + Percentage_forest, data = umf, keyfun = "halfnorm", output = "abund") Abundance (log-scale): Estimate SE z P(>|z|) (Intercept) -1.225 0.229 -5.34 9.47e-08 April_temp_degrees_celsius -0.615 0.139 -4.43 9.33e-06 Percentage_forest 0.574 0.112 5.12 3.10e-07 Detection (log-scale): Estimate SE z P(>|z|) 3.28 0.138 23.8 2.88e-125 AIC: 371.3509 Number of sites: 201 Survey design: line-transect Detection function: halfnorm UnitsIn: m UnitsOut: ha

Thanks, 

Christina

[Ken Kellner]

Hi Christina,

This was also unintuitive to me, but after some simulations (attached), I think your result is expected.

I think what is happening here, crudely speaking, is that the detection submodel drives the overall abundance, and then within that constraint the site-level abundance/density covariates control the relative abundance prediction at each site. Someone else on the list might be able to give a better explanation.

As you may have noticed, if you change the covariates on sigma, the predict sum does change.

For getting an estimate of total abundance across sites, it is likely better to use the ranef() and bup() functions instead of predict. These functions give an estimate of the latent abundance at each site, incorporating information about what was actually observed at each site. These values will not necessarily be the same for different models based on the simulation.

Ken

> --
> *** Three hierarchical modeling email lists ***
> (1) unmarked (this list): for questions specific to the R package unmarked
> (2) SCR: for design and Bayesian or non-bayesian analysis of spatial capture-recapture
> (3) HMecology: for everything else, especially material covered in the books by Royle & Dorazio (2008), Kéry & Schaub (2012), Kéry & Royle (2016, 2021) and Schaub & Kéry (2022)
> ---
> You received this message because you are subscribed to the Google Groups "unmarked" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to unmarked+u...@googlegroups.com.
> To view this discussion visit https://groups.google.com/d/msgid/unmarked/545bf2ee-a4f5-489c-b645-7ca76a6cb25bn%40googlegroups.com.

[Christina Hunt]

Hi Ken,

Thanks very much for looking into this for me and for producing that simulation code. It is good to know that my results are to be expected and are nothing to worry about. 

Yes I had noticed that the predict sum changes when I change the covariates on sigma. It seems a bit surprising to me that the overall abundance is driven only by the detection model, rather than a combination of detection and abundance, but at least now I know that is what's expected with these models and isn't something that I have done wrong.

Thank you as well for the advice on using bup() and ranef() rather than predict() for my abundance estimates.

Thanks,

Christina