I'm modeling counts of organisms over time at eleven locations. I'd like to account for temporal autocorrelation in the counts, assuming it's present.
As the data are not equally spaced in time, I'm exploring empirical variograms (and comparing variogram model fits) based on the residuals of a generalized linear model assuming independence.
As shown below, the variogram based on the classical semivariance estimator looks quite different from the robust (Cressie) estimator. Furthermore, a fit of most variogram models to the classical variogram produces singularity errors (with good reason based on the variogram plot) while a decent exponential model (e.g.) can be fitted to the robust variogram. Unfortunately, it seems as though
nlme:::lme (which I'm using for the generalized model) only uses the classical variogram when fitting correlation structures.
Thus, my questions:
(1) Is the use of the robust estimator for constructing the correlation structure justified in this case?
(2) If so, is it appropriate to fit the the robust variogram model in, e.g.,
geoR, and then specify/fix the correlation structure in the
# Download data datURL <- "https://www.dropbox.com/s/54q9ocwztt3swap/rsiggeo_variogram.csv" dat <- repmis::source_data(datURL, sep = ",", header = TRUE) # Load necessary packages library(gstat) library(sp) # Convert to appropriate class # Y coordinates are arbitrarily assigned values for each of 11 sites # (count locations; assumed independent); Y coords were assigned such # that each site is separated from others by more than the cutoff distance # used in variogram construction coordinates(dat) <- ~x+y # Classic variogram model with example exponential fit varClassic <- variogram(resid ~ 1, dat, cutoff = 180, width = 180/25) # Note singularity error; starting values determine 'fit' expClassic <- fit.variogram(varClassic, model = vgm(1, "Exp", 1, 0)) plot(varClassic, expClassic, xlab = 'Distance (days)')
# Robust (Cressie) variogram model with example exponential fit varCressie <- variogram(resid ~ 1, dat, cutoff = 180, width = 180/25, cressie=TRUE) expCressie <- fit.variogram(varCressie, model = vgm(0.3, "Exp", 50, 0.1)) plot(varCressie, expCressie, xlab = 'Distance (days)')