Question about making predictions when using a hurdle approach

[sashavas3]

Hi,

I'm using a hurdle approach to assess canopy insect responses to some environmental variables. The insects are classified as generalists, specialists and omnivores (the trait included in the models). I am wondering though if I'm using the correct approach to make the predictions.

So far I've run the codes below to obtain predictions for total richness (from the probit model) and total abundance (from the abundance conditional on presence model), and for the community weighted mean trait values for each model. Is it correct to obtain separate predictions for richness and abundance conditional on presence?

I ask because in the HMSC book there is a section that states:
"Once the two components of the hurdle model have been fitted, they can be combined to produce
predictions that predict the full distribution of counts in the original response variable, including the
zeros and non-zeros. To do so, first the presence-absence model is used to predict whether the
species is present or not. If it is predicted to be absent, then the prediction is zero, and the model for
abundance (conditional on presence) is not needed at all. If the species is predicted to be present,
then the model for abundance (conditional on presence) is used to predict how many individuals
there are."

I haven't been able to find an example of how to do this in the book or in the vignettes, and I've also tried to see how it has been done in scripts provided in other papers, but so far haven't found an example.

Is the approach I'm using correct, or should I definitely combine the models to make predictions? And if so, how can I do this?

Any advice or suggestions are very welcome!

Sasha

Total response to temperature

Gradient = constructGradient(hM = model, focalVariable ="Temperature", non.focalVariables = 1)
predY = predict(model, XData = Gradient$XDataNew, studyDesign = Gradient$studyDesignNew, ranLevels = Gradient$rLNew, expected = TRUE)
png(sprintf("canopy-temp.png"), width=11.5, height=12, unit ="cm", res=600)
plotGradient(model, Gradient, pred = predY, measure = "S", index = 2, showData = TRUE,
ylabel = "Total richness", xlabel = "Temperature",
cicol = rgb(0.788,0.788,0.788, alpha = 0.5), showPosteriorSupport = FALSE,
q = c(0.08,0.5,0.92)) # 84% CI
dev.off()

Omnivore response to temperature

Gradient = constructGradient(hM = model, focalVariable ="Temperature", non.focalVariables = 1)
predY = predict(model, XData = Gradient$XDataNew, studyDesign = Gradient$studyDesignNew, ranLevels = Gradient$rLNew, expected = TRUE)
png(sprintf("canopy-omnivore-temp.png"), width=11.5, height=12, unit ="cm", res=600)
plotGradient(model, Gradient, pred = predY, measure = "T", index = 2, showData = TRUE,
ylabel = "Omnivores", xlabel = "Temperature",
cicol = rgb(0.788,0.788,0.788, alpha = 0.5), showPosteriorSupport = FALSE,
q = c(0.08,0.5,0.92)) # 84% CI
dev.off()

Specialist response to temperature

Gradient = constructGradient(hM = model, focalVariable ="Temperature", non.focalVariables = 1)
predY = predict(model, XData = Gradient$XDataNew, studyDesign = Gradient$studyDesignNew, ranLevels = Gradient$rLNew, expected = TRUE)
png(sprintf("canopy-specialist-temp.png"), width=11.5, height=12, unit ="cm", res=600)
plotGradient(model, Gradient, pred = predY, measure = "T", index = 3, showData = TRUE,
ylabel = "Specialists", xlabel = "Temperature",
cicol = rgb(0.788,0.788,0.788, alpha = 0.5), showPosteriorSupport = FALSE,
q = c(0.08,0.5,0.92)) # 84% CI
dev.off()**