How to deal with spatially autocorrelated residuals in GLMM

user_10116914 0
Viewed 1782

I am conducting an analysis of where on the landscape a predator encounters potential prey. My response data is binary with an Encounter location = 1 and a Random location = 0 and my independent variables are continuous but have been rescaled.

I originally used a GLM structure 

glm_global <- glm(Encounter ~ Dist_water_cs+coverMN_cs+I(coverMN_cs^2)+
                  Prey_bio_stand_cs+Prey_freq_stand_cs+Dist_centre_cs,
                  data=Data_scaled, family=binomial)

but realized that this failed to account for potential spatial-autocorrelation in the data (a spline correlogram showed high residual correlation up to ~1000m). 

Correlog_glm_global <- spline.correlog (x = Data_scaled[, "Y"],
                                        y = Data_scaled[, "X"],
                                        z = residuals(glm_global,
                                        type = "pearson"), xmax = 1000)

I attempted to account for this by implementing a GLMM (in lme4) with the predator group as the random effect. 

glmm_global <- glmer(Encounter ~ Dist_water_cs+coverMN_cs+I(coverMN_cs^2)+
                     Prey_bio_stand_cs+Prey_freq_stand_cs+Dist_centre_cs+(1|Group),
                     data=Data_scaled, family=binomial)

When comparing AIC of the global GLMM (1144.7) to the global GLM (1149.2) I get a Delta AIC value >2 which suggests that the GLMM fits the data better. However I am still getting essentially the same correlation in the residuals, as shown on the spline correlogram for the GLMM model).

Correlog_glmm_global <- spline.correlog (x = Data_scaled[, "Y"],
                                         y = Data_scaled[, "X"],
                                         z = residuals(glmm_global,
                                         type = "pearson"), xmax = 10000)

I also tried explicitly including the Lat*Long of all the locations as an independent variable but results are the same. After reading up on options, I tried running Generalized Estimating Equations (GEEs) in “geepack” thinking this would allow me more flexibility with regards to explicitly defining the correlation structure (as in GLS models for normally distributed response data) instead of being limited to compound symmetry (which is what we get with GLMM). However I realized that my data still demanded the use of compound symmetry (or “exchangeable” in geepack) since I didn’t have temporal sequence in the data. When I ran the global model

gee_global <- geeglm(Encounter ~ Dist_water_cs+coverMN_cs+I(coverMN_cs^2)+
                     Prey_bio_stand_cs+Prey_freq_stand_cs+Dist_centre_cs,
                     id=Pride, corstr="exchangeable", data=Data_scaled, family=binomial)

(using scaled or unscaled data made no difference so this is with scaled data for consistency) suddenly none of my covariates were significant. However, being a novice with GEE modelling I don’t know a) if this is a valid approach for this data or b) whether this has even accounted for the residual autocorrelation that has been evident throughout.

I would be most appreciative for some constructive feedback as to 1) which direction to go once I realized that the GLMM model (with predator group as a random effect) still showed spatially autocorrelated Pearson residuals (up to ~1000m), 2) if indeed GEE models make sense at this point and 3) if I have missed something in my GEE modelling. Many thanks.

geospatial data-modeling r
1 Answers

Taking the spatial autocorrelation into account in your model can be done is many ways. I will restrain my response to R main packages that deal with random effects.

First, you could go with the package nlme, and specify a correlation structure in your residuals (many are available : corGaus, corLin, CorSpher ...). You should try many of them and keep the best model. In this case the spatial autocorrelation in considered as continous and could be approximated by a global function.

Second, you could go with the package mgcv, and add a bivariate spline (spatial coordinates) to your model. This way, you could capture a spatial pattern and even map it. In a strict sens, this method doesn't take into account the spatial autocorrelation, but it may solve the problem. If the space is discret in your case, you could go with a random markov field smooth. This website is very helpfull to find some examples : https://www.fromthebottomoftheheap.net

Third, you could go with the package brms. This allows you to specify very complex models with other correlation structure in your residuals (CAR and SAR). The package use a bayesian approach.

I hope this help. Good luck

user_53088210
Related