When applying their principal coordinates of neighbor matrices (PCNM; Borcard et al. 2002, Ecological Modelling) method to real ecological data, Borcard et al. (2004, Ecology) recommended to check for linear spatial trends in the response data (linear gradients), prior to the analysis. They advocated that linear trends in the response data could obscure other recoverable structures in the data at finer scales.

When Dray et al. (2006, Ecological Modelling) introduced their Moran's Eigenvector Maps (MEMs), they showed that PCNMs are a particular case of MEMs. When they analysed one of the datasets of Borcard et al. (2004; the "famous" oribatid mites dataset), they also detrended the response variables (species abundances) by multiple linear regression on geographic coordinates. Using this same dataset, Jombart et al. (2009, Ecography) also detrended the environmental variables.

However, in my dataset (irregular 2D sampling grid), regression on geographic coordinates both: 1) prevents from selecting the best method for constructing the spatial weighting matrix, and 2) when the regression line (e.g. abundance vs. coordinate x) has a negative intercept, the regression can have positive residuals at abundance values = 0, for small-enough values of the geographic coordinates. This hampers the interpretation of the analysis.

Is it mandatory to detrend both the response variables and the environmental variables, if a linear gradient is found? In this case, how to solve the mentioned issues?

We're working in R, adespatial package (https://cran.r-project.org/web/packages/adespatial/vignettes/tutorial.html#swm)