predict() when using a pre-scaled variable and quadratic term
Tim Doherty
Hi
I am getting different predictions for two models where the model summaries are exactly the same. The models include a quadratic term. It must be to do with the way I am specifying the newdata or using predict(), but I can't figure it out. Any help is greatly appreciated!
Important background is that I am scaling my predictor variable (YSF) prior to fitting the model because when I used scale() within the model call, predictions did not work when for second variable (not included in the models below). Something similar was described in an older post on this forum.
With regards to the code below, all of the values from summary() are exactly the same for occ1.1 and occ1.2, but when I predict using new data, the predictions are very different. I note that this does not occur if I remove the quadratic term from the model and just fit YSF.
Can you see anything wrong with the way I am specifying the newdata or using predict()?
I have checked that the scaled values in newdata1.2 match with the unscaled values for newdata1.1. I also tried changing the order of scaling and transforming, but that didn't help.
Many thanks.
# scale YSF for model 1.2
umf_list$Woyl@siteCovs$YSF_sc = scale(umf_list$Woyl@siteCovs$YSF)
# store the YSF transformation values
ysf_center <- attr(scale(umf_list$Woyl@siteCovs$YSF), "scaled:center")
ysf_scale <- attr(scale(umf_list$Woyl@siteCovs$YSF), "scaled:scale")
# fit model 1.1 where scaling is done inside the model
occ1.1 = occu(~Track + Season
~scale(YSF) + I(scale(YSF)^2) +
(1|LocationName_numeric) + (1|Year),
data = umf_list$Woyl)
# fit model 1.2 with the pre-scaled YSF
occ1.2 = occu(~Track + Season
~YSF_sc + I(YSF_sc^2) +
(1|LocationName_numeric) + (1|Year),
data = umf_list$Woyl)
# predict the first model (occ1.1)
newdata1.1 = data.frame(YSF = 0:55)
pred1.1 = unmarked::predict(occ1.1, newdata=newdata1.1,
type="state", append=TRUE, re.form = NA)
# predict the second model (occ1.2)
YSF_seq = 0:55
newdata1.2 <- data.frame(YSF_sc = (YSF_seq - ysf_center) / ysf_scale)
pred1.2 = predict(occ1.2, newdata=newdata1.2,
type="state", append=TRUE, re.form = NA)
Ken Kellner
Nesting multiple functions inside formulas, e.g. scale(I(x^2)) or I(scale(x)^2), doesn't work properly in unmarked with predict. The models give the correct results (which is why your two models have identical output) but when you try to predict with the nested formula model you get the wrong answer because the scaling is not re-applied correctly. So, I would use the pre-scaled model. I would recommend calculating a few predictions manually to make sure, though.
The last time I looked into this I determined the issue was happening in the underlying model.matrix function in R, which meant it isn't something we could fix in unmarked. Probably we should detect this nested pattern in the formula and issue a warning.
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/7d5cf992-f592-463f-b69f-a6f990d532b4n%40googlegroups.com.