### copy data into 'dat' and examine data
dat <- dat.maire2019$dat
dat[-10]
### copy distance matrix into 'dmat' and examine first 5 rows/columns
dmat <- dat.maire2019$dmat
dmat[1:5,1:5]
if (FALSE) {
### load metafor package
library(metafor)
### fit a standard random-effects model ignoring spatial correlation
res1 <- rma.mv(s1, vars1, random = ~ 1 | site_station, data=dat)
res1
### fit model allowing for spatial correlation
res2 <- rma.mv(s1, vars1, random = ~ site_station | const, struct="SPGAU",
data=dat, dist=list(dmat), control=list(rho.init=10))
res2
### add random effects for sites and stations within sites
res3 <- rma.mv(s1, vars1, random = list(~ 1 | site/station, ~ site_station | const), struct="SPGAU",
data=dat, dist=list(dmat), control=list(rho.init=10))
res3
### likelihood ratio tests comparing the models
anova(res1, res2)
anova(res2, res3)
### profile likelihood plots for model res2
profile(res2, cline=TRUE)
### effective range (river-km for which the spatial correlation is >= 0.05)
sqrt(3) * res2$rho
### note: it was necessary to adjust the starting value for rho in models
### res2 and res3 so that the optimizer does not get stuck in a local maximum
profile(res2, rho=1, xlim=c(0,200), steps=100)
}
Run the code above in your browser using DataLab