I use gstat to predict a binomial data, but the predicted values go above 1 and below 0. Does anyone know how I can deal with this issue? Thanks.
data(meuse)
data(meuse.grid)
coordinates(meuse) <- ~x+y
coordinates(meuse.grid) <- ~x+y
gridded(meuse.grid) <- TRUE
#glm model
glm.lime <- glm(lime~dist+ffreq, meuse, family=binomial(link="logit"))
summary(glm.lime)
#variogram of residuals
var <- variogram(lime~dist+ffreq, data=meuse)
fit.var <- fit.variogram(var, vgm(nugget=0.9, "Sph", range=sqrt(diff(meuse@bbox\[1,\])^2 + diff(meuse@bbox\[2,\])^2)/4, psill=var(glm.lime$residuals)))
plot(var, fit.var, plot.nu=T)
#universal kriging
kri <- krige(lime~dist+ffreq, meuse, meuse.grid, fit.var)
spplot(kri[1])
In general, with this kind of regression kriging approach there is no guarantee that the model will be valid as the calculation of the trend and the residuals is separated. A few notes on your code. Notice that you use
variogram
to calculate the residual variogram, butvariogram
uses a normal linear model to calculate the trend and thus also the residuals. You need to determine your residuals from yourglm
, and then calculate a residual variogram based on that.You could do this manually, or have a look at the
fit.gstatModel
function from theGSIF
package. You could also have a look atbinom.krige
from thegeoRglm
package. This thread on R-sig-geo might also be interesting: