A fellow R-sig-eco list member assisted me in working out a method for dealing with spatial autocorrelation in classification trees (conditional inference trees in R-package party).

The whole idea is actually plain simple: compute the classification tree, calculate residuals and use it for a Mantel-test and Mantel correlograms.

The Mantel correlograms test differrences in dissimilarities of

the residuals across several spatial distances and thus enable you to detect lag-distances where possible spatial autocorrelation vanishes.

I did the mantel tests and correlograms with the package ecodist.

The below toy example is obviously not spatially autocorrelated - as is shown by the Mantel tests and correlograms. If I'd encounter autocorrelation with real-world data, I'd try to use subsamples of the data avoiding resampling within the lag-distance..

For an example with this scheme in action see this paper by Bálint Czúcz et al.

library(party) # for simplicity only one env-factor: env <- data.frame(F1 = as.factor(rep(LETTERS[1:4], each = 40))) n <- nrow(env) # probs for simulated sampling: p <- seq(0, 1, 1/n) # simulated binary response: binresp <- numeric() for(i in 1:n) { binresp[i] <- sample(0:1, 1, replace = T, c(p[i], 1-p[i])) } binresp <- as.factor(binresp) ct <- ctree(binresp ~ env$F1) plot(ct) # predicted nodes: pred <- names(predict(ct)) # predicted probabilties: ct_resp <- treeresponse(ct) # pull second elements of list items: ct_resp <- sapply(ct_resp, "[[", 2) # convert response from factor back to numeric: binresp <- as.numeric(ifelse(binresp == 1, 1, 0)) # calculate residuals: Residuals <- ct_resp - binresp cbind(ct_resp, binresp, Residuals) # simulate XY-coords: Y <- jitter(sample(39:42, n, replace = T)) X <- jitter(sample(10:12, n, replace = T)) XY <- cbind(X, Y) # install.packages("ecodist") library(ecodist) # make distance matrices: d_XY <- lower(dist(XY)) d_Resid <- lower(dist(Residuals)) # check distribution of residuals: table(d_Resid); hist(d_Resid) # mantel tests (with Pearson's and Spearman's Corr. Coeff.): mantel(d_Resid ~ d_XY, nperm = 10000) mantel(d_Resid ~ d_XY, nperm = 10000, mrank = T) # mantel correlograms: plot(mgram(d_Resid, d_XY)) plot(mgram(d_Resid, d_XY, mrank = T))

## No comments :

## Post a Comment