continuous covariate on density and predictions

[Gilles Maurer]

Dear all,

I am struggling to get density predictions using continuous covariate on density, as well as getting the total N. 
It worked with a simplified statespace (at low resolution) for testing purposes.  
Since, I build a statespace including non habitat and a huge number of cell grids : 
- sspace : around 100k cells, each cell with its cov values,  like scaled elevation (elev) in the sspace.

- sf :  2386 traps with capthistory. These traps are indeed grid cells visited for feces collection.  elev is also part of the trapCovs
- one session & 7 occasions

Null model is fine (D ~ 1 p0 ~ 1 sig ~ 1)
I got real estimates for density, detection, sigma, total N, realzed density etc..

Model with elev as a covariate of density is also running (D ~ occupancy p0 ~ 1 sig ~ 1) but takes some time to compute using trimS ! 
Getting real estimates for detection and sigma is OK 
But not convenient for density as I ended with a table of 200k rows ( 100k per sex) with an estimate, se, etc for each grid cell per sex.

pred.df.dens <- data.frame(Session = factor(1),
                           elev = s_space[[1]]$elev)

#make predictions on the real scale
pred.dens <- get.real(model, type = "dens",
                      newdata = pred.df.dens,
                      d.factor = d.factor)

# make predictions on the real scale - total N
(pred.dens <-get.real(model,type ="dens",newdata =pred.df.dens,d.factor =nrow(model$ssDF[[1]])))

How can I get the total N easily ? How to run AC per indivual ? Realized density ?
I guess I can plot the pred.dens as a raster for  density in sqkm . 

 
PS : I did not find the script of the ocelot or Bear example to help

Thank you for your help
Gilles

[Daniel Linden]

Hi Gilles, I'm confused when you say that the null model allowed you to get total N and realized density.  How did you calculate those?

The predict.oSCR function should provide what you need, either way.

[Gilles Maurer]

Thanks Dan for your answer

I realized that maybe I was not clear between get.real function on the model results and predict.oSCR ! 

What I did on  the following model  : 

model0wide <- oSCR.fit(model = #model formulation
                     list(D ~ 1,   p0 ~ 1,      sig ~ 1),     scrFrame = sf,   ssDF = data.ss, trimS=10000)

1 - make predictions on the real scale  
 
pred.df.dens <-data.frame(Session =factor(1))

#make predictions on the real scale

(pred.dens <-get.real(model0wide,type ="dens",newdata =pred.df.dens, d.factor=d.factor))

and 

2 - the predict function

pred <- predict.oSCR(model0wide, sf, data.ss, override.trim=TRUE) 

And it works well,

Then I added a covariate to the density on a new 'model'  oSCR.fit(model = list(D ~ elev, p0 ~ 1, sig ~ 1), scrFrame = sf, ssDF = s_space, trimS = 15000)

I can run the predict.oSCR function and get for instance the prob for each  indiv accross the study site (from  pred$preds[[1]] ) 
I also run the get.real function on sigma or detection, it works well as I do not have covariates

For density : 

I have one value of elevation (elev) for each cell of the statespace.  

# Build the prediction matrix accross statespace cells 

pred.df.dens <- data.frame(Session = factor(1), elev = s_space[[1]]$elev)

=> the dataframe has twice more lines ( male + female) than the number of cells in the statespace. = 2 x 100 000 rows

   
#make predictions on the real scale
pred.dens <- get.real(model, type = "dens",
                      newdata = pred.df.dens,
                      d.factor = d.factor)

I obtain a table of 200 000 rows with a density per sex per pixel while when using the same code with the previous model I got total of male and total of female.

Hope this is more clear, 

Have a nice week end 

Gilles

[Daniel Linden]

Hi Gilles,

The get.real function will provide a density prediction for each pixel in your newdata. Since you asked get.real to use a data frame with a single record (indicating Session 1) for your newdata, your prediction was for 2 values (male and female expected densities).  The d.factor is a scalar that gets multiplied by the pixel prediction so when d.factor is >1, the predictions no longer represent a state space pixel prediction.  In your case, if d.factor was the number of pixels in the state space, your estimates now represent the total population sizes in the state space.

In your D~1 model, if you had used get.real with your original state space as the newdata, you would get an estimate for every pixel.  The sum of those pixel density estimates would give you exactly the same value as the single-record data frame you provided.  So, for your D~elev model, you can use the following to estimate total population sizes:

pred.dens <- get.real(model, type = "dens",
                      newdata = pred.df.dens,

                      d.factor = 1)

Note, d.factor is 1 by default so you don't have to specify but I'm emphasizing that it needs to be 1 for this to work.  You then sum the estimates for each sex to generate population size.  Error around those estimates is more difficult, but Chris has been working on a function for that.

Let me know if that works!

[Gilles Maurer]

Thanks Dan for your answer,

it worked perfectly, but  cannot get error on abundance by suming density per pixel. I know that SCR is mainly focusing on density rather than abundance. Density per pixel and associated error is fine for me, but conservationists would also like to get abundance estimates. 
  You mentioned that Chris has been working on a function for that purpose. Any feedback on that ? 
Is bootstrapping an option ? Boostrapping on d.beta would be more robust than bootstrapping prediction/pixel, but it seems very tricky
Any thoughts ?

Have a nice day and thanks again for your help
Gilles 
 

[Daniel Linden]

Hi Gilles, I'll let Chris chime in on his function, he can probably share it directly with you before we add it to the package.

Bootstrapping is always a possibility, but the main reason we have not used it with oSCR is the prohibitive computational demand (depending on your data dimensions and model complexity).  I have typically transitioned to MCMC with something like NIMBLE for this type of task because if I'm going to wait for bootstrap results, I would rather just have posterior distributions.  Obviously, it's not always that straightforward...

[Gilles Maurer]

Dear Dan, 
thanks again for your answer, 

If the function could be available sometimes, I will be happy to try, otherwise MCMC could be a great option, 

have a nice day

Gilles

--
You received this message because you are subscribed to a topic in the Google Groups "oSCR" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/oscr_package/TEG7O8eD_aU/unsubscribe.
To unsubscribe from this group and all its topics, send an email to oscr_package...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/oscr_package/0ef13085-ea7b-4737-a48e-cb44674cc12en%40googlegroups.com.