I have an ecological data set, whereby sediment mud content (%) (i.e. the continuous explanatory variable) is thought to be explaining the spatial distribution (i.e. presence/absence) of various estuarine benthic macroinvertebrate taxa (i.e. the binary response variable). Below is a subset of the data (actual data consists of >1400 replicates) for a single taxon, where the response is 'Bin' and the single predictor is 'Mud' and 'Estuary' is considered a random effect.
Estuary Mud Bin FW 0.1 0 FW 0.3 0 FW 0.2 0 FW 2.1 0 POR 7.1 0 FOR 11.4 0 JRE 25.5 1 JRE 33.6 1 JRE 44.6 1 JRE 55.6 1 JRE 63.6 1 FW 76.6 1 FW 87.6 1 FW 90.6 1 FW 93.6 1
Here's the GAM model which used the
mgcv package (and was conditioned on the 1400+ data points):
aa1.estuary<-gam(Bin~s(Mud,bs="ps", k=6) + s(Estuary, bs="re"),family=binomial, gamma=1,data=sample)
Here's the Cross validation step, which used
CVgam in the
gamclass package, but returned the error below when
Estuary was included as a random effect:
aa1.val <- CVgam(formula=Bin~s(Mud, bs="ps", k=6) + s(Estuary, bs="re"), data = sample, nfold = 10, debug.level = 0, printit = TRUE, method = "GCV.Cp",cvparts = NULL, gamma = 1, seed =100) Error in X[ind, ] : subscript out of bounds
Estuary was omitted, generating the following output:
aa1.val <- CVgam(formula=Bin~s(Mud, bs="ps", k=6) , data = sample, nfold = 10, debug.level = 0, printit = TRUE, method = "GCV.Cp",cvparts = NULL, gamma = 1, seed =100) GAMscale CV-mse-GAM 0.1491 0.1499
Hence I have two related questions: (1) whats causing the error when
Estuary is included in the cross-validation step? and (2) do the above estimates of
CV-mse-GAM =0.1499seem too close to one another to be realistic?
Thanks in advance,