Accounting for spatial autocorrelation and grouping random effect in nlme

Question

I would like to test differences in mean values of a Climate_variable between 4 different climatic Scenarios that are projected for many (3000) surveyed Sites for which I know the coordinates (Lon and Lat). My database is structured as follows:

Climate_variable	Scenario	Site	Lon	Lat
0.85	A	S1	19.37	41.85
1.13	B	S1	19.37	41.85
1.55	C	S1	19.37	41.85
2.12	D	S1	19.37	41.85
0.74	A	S2	21.42	40.36
1.08	B	S2	21.42	40.36
1.44	C	S2	21.42	40.36
2.01	D	S2	21.42	40.36

Because I have 4 observations for each Site but I am not interested in this effect, I wanted to go for a Linear Mixed Model with Site as random effect. However, climatic variables are often highly spatially autocorrelated so I also wanted to add a spatial autocorrelation structure using the coordinates of the sites.

This is the code I run using the nlme package on R:

mod <- lme(Climate_variable ~ Scenario,
           correlation = corExp(form = ~Lon+Lat),
           random = list(~1|Site))

However, there is an error message saying cannot have zero distances in "corSpatial", I guess because I have 4 observation with the same coordinates for each site.

So I tried this instead:

mod <- lme(Climate_variable ~ Scenario,
           correlation = corExp(form = ~Lon+Lat | Scenario),
           random = list(~1|Site))

But I get another error message: incompatible formulas for groups in 'random' and 'correlation'

How could I account for both the spatial autocorrelation and the non-independance of the observations for different Scenarios on the same Site?

Answer 1

I have to tell you that i do it different. First i make the null model indicating fixed and random. If you put site as random, you are telling the program to control for the different between sites, that i think is not what you want. right? You want to control for the difference between the scenarios. Then scenario is random. example from the package.

null.model <- lme(fixed = thick ~ 1, data = spdata, random = ~ 1 | dummy)

exp.sp <- update(null.model, correlation = corExp(1, form = ~ east + north), method = "ML")

Your error."There is an error message saying cannot have zero distances in "corSpatial", I guess because I have 4 observation with the same coordinates for each site." Exactly is because you have the same coordinates and when the program try to calculate the geographic distance to check for autocorrelation give 0, that for the program sound absurd. I think this method is not the correct if you only have one site. This is for checking spatial autocorrelation. It is better you only check for correlation in your scenarios in a simple way.

Answer 2

Even though this analysis might not be relevant (see comment from @mkt), this code works:

mod <- lme(Climate_variable ~ Scenario, correlation = corExp(form = ~(Lon+Lat)/Scenario), random = list(~1|Site))