R/msr.mantelrtest.R

In adespatial: Multivariate Multiscale Spatial Analysis

#' Moran spectral randomization for Mantel test
#'
#' This function allows to test the Mantel statistic using constrained null
#' models in the presence of spatial autocorrelation. Random replicates of the
#' second distance matrix are produced. They are spatially-constrained to
#' preserve the global autocorrelation (Moran's I) and the spatial structures at
#' multiple scales. Multiscale property is defined by the power spectrum (i.e.
#' decomposition of the variance of the original variables) on a basis of
#' orthonormal eigenvectors (Moran's Eigenvector Maps, MEM).
#'
#' @param x An object generated by the \code{mantel.randtest} function.
#' @param listwORorthobasis an object of the class \code{listw} (spatial
#'   weights) created by the functions of the \pkg{spdep} package or an object
#'   of class \code{orthobasis}
#' @param nrepet an \code{integer} indicating the number of replicates
#' @param method an character specifying which algorithm should be used to
#'   produce spatial replicates (see \code{\link{msr.default}}).
#' @param \dots further arguments of the \code{\link{msr.default}} function.
#' @returns An object of class \code{randtest}.
#' @author Sylvie Clappe, Stephane Dray \email{stephane.dray@@univ-lyon1.fr}
#' @seealso \code{\link{msr.default}}, \code{\link[ade4]{mantel.randtest}}
#' @references Crabot, J., Clappe, S., Dray, S. and Datry, T. (2019) Testing the
#'   Mantel statistic with a spatially-constrained permutation procedure.
#'   Methods in Ecology and Evolution. in press.
#'
#' @keywords spatial
#' @examples
#'
#' if(require("ade4", quietly = TRUE)
#' & require("spdep", quietly = TRUE)){
#' data(mafragh, package = "ade4")
#'
#' d1 <- dist(mafragh$env[,1:3])
#' d2 <- dist(mafragh$env[,7])
#' t1 <- mantel.randtest(d1,d2)
#' t1
#'
#' lw <- nb2listw(mafragh$nb)
#' t2 <- msr(t1, listwORorthobasis = lw)
#' t2
#'
#' }
#'
#' @importFrom ade4 pcoscaled is.euclid dist.quant
#' @export
msr.mantelrtest <- function(x, listwORorthobasis, nrepet = x$rep, method = c("pair", "triplet", "singleton"), ...){
    
    appel <- as.list(x$call)
    dX <- eval.parent(appel$m1)
    dY <- eval.parent(appel$m2)
    if(!is.euclid(dX)){
        stop("This test only works for Euclidean distances (m1 is not).")
    }
    X <- pcoscaled(dX)
    obs <- cor(dist(X), dY)
    Xmsr <- msr(X, listwORorthobasis, nrepet = nrepet, method = method, simplify = FALSE, ...)
    
    sim <- sapply(Xmsr, function(i) cor(dist(i), dY)) 
    
    test <- as.randtest(obs = obs, sim = sim, 
        call = match.call())
    return(test)
}