Moran spectral randomization

Source: R/msr.R

This function allows to generate spatially-constrained random variables preserving 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). The function provides methods for univariate randomization, joint randomization of a group of variables while keeping within-group correlations fixed and univariate randomization with a fixed correlation between original data and randomized replicates.

msr(x, ...)

# Default S3 method
msr(
  x,
  listwORorthobasis,
  nrepet = 99,
  method = c("pair", "triplet", "singleton"),
  cor.fixed,
  nmax = 100,
  simplify = TRUE,
  ...
)

Arguments

x

For msr.default, a vector, a matrix or a data.frame with the original variables. If NCOL(x) > 1, then the joint randomization procedure that preserves the correlations among variables is used.

...

further arguments passed to or from other methods

listwORorthobasis

an object of the class listw (spatial weights) created by the functions of the spdep package or an object of class orthobasis

nrepet

an integer indicating the number of replicates

method

an character specifying which algorithm should be used to produce spatial replicates (see Details).

cor.fixed

if not missing, the level of correlation between the original variable and its randomized replicates

nmax

the number of trials used in the "triplet" procedure.

simplify

A logical value. If TRUE, the outputs for univariate procedures are returned in a matrix where each column corresponds to a replicate. If FALSE a list is returned.

Value

Either a matrix (if simplify is TRUE) or a list with randomized replicates.

Details

Three procedures are implemented in the function. The "pair" procedure is the more general as it can be applied in the three cases (univariate, univariate with fixed correlation and multivariate). This procedure preserves the power spectrum by pair of MEMs but not strictly the global autocorrelation level (Moran's I). The "singleton" procedure can be used for univariate and multivariate cases. It preserves strictly the global level of autocorrelation and the power spectrum. The "triplet" procedure can only be applied in the univariate case. It preserves the power spectrum by triplet of MEMs and strictly the global autocorrelation level.

References

Wagner, H.H. and Dray S. (2015) Generating spatially-constrained null models for irregularly spaced data using Moran spectral randomization methods. Methods in Ecology and Evolution, 6: 1169-1178. doi:10.1111/2041-210X.12407

Author

Stephane Dray stephane.dray@univ-lyon1.fr and Helene H Wagner helene.wagner@utoronto.ca

Examples


library(spdep)
x1 <- matrix(rnorm(81*5), nrow = 81)
lw1 <- nb2listw(cell2nb(9, 9))

moran.mc(x1[,1], lw1, 2)$statistic
#>  statistic 
#> 0.08901257 

## singleton
x1.1 <- msr(x1[,1], lw1, nrepet = 9, method = "singleton")
apply(x1.1, 2, function(x) moran.mc(x, listw = lw1, nsim = 2)$statistic)
#> [1] 0.08901257 0.08901257 0.08901257 0.08901257 0.08901257 0.08901257 0.08901257
#> [8] 0.08901257 0.08901257

## triplet
x1.2 <- msr(x1[,1], lw1, nrepet = 9, method = "triplet")
apply(x1.2, 2, function(x) moran.mc(x, listw = lw1, nsim = 2)$statistic)
#> [1] 0.08901257 0.08901257 0.08901257 0.08901257 0.08901257 0.08901257 0.08901257
#> [8] 0.08901257 0.08901257

## pair
x1.3 <- msr(x1[,1], lw1, nrepet = 9, method = "pair")
apply(x1.3, 2, function(x) moran.mc(x, listw = lw1, nsim = 2)$statistic)
#> [1] 0.08853303 0.09114309 0.09651904 0.08944276 0.08881118 0.09331037 0.08668055
#> [8] 0.09293674 0.09300961

## pair with cor.fixed
x1.4 <- msr(x1[,1], lw1, nrepet = 9, cor.fixed = 0.5)
apply(x1.4, 2, function(x) moran.mc(x, listw = lw1, nsim = 2)$statistic)
#> [1] 0.09022067 0.09567686 0.08808840 0.09551797 0.08840153 0.09349472 0.09440850
#> [8] 0.08910824 0.08807222
cor(x1[,1], x1.4)
#>      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
#> [1,]  0.5  0.5  0.5  0.5  0.5  0.5  0.5  0.5  0.5

## pair preserving correlations for multivariate data
x1.5 <- msr(x1, lw1, nrepet = 9, cor.fixed = 0.5)
cor(x1)
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
lapply(x1.5, cor)
#> [[1]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[2]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[3]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[4]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[5]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[6]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[7]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[8]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[9]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 

apply(x1, 2, function(x) moran.mc(x, listw = lw1, nsim = 2)$statistic)
#> [1]  0.08901257 -0.01819323  0.03040565 -0.02081023  0.01448972
apply(x1.5[[1]], 2, function(x) moran.mc(x, listw = lw1, nsim = 2)$statistic)
#> [1]  0.08883646 -0.02639149  0.02498800 -0.02419777  0.01451139

## singleton preserving correlations for multivariate data
x1.6 <- msr(x1, lw1, nrepet = 9, method = "singleton")
cor(x1)
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
lapply(x1.6, cor)
#> [[1]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[2]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[3]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[4]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[5]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[6]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[7]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[8]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 
#> [[9]]
#>             [,1]        [,2]        [,3]        [,4]        [,5]
#> [1,]  1.00000000 -0.11247063 -0.02124780  0.06779336 -0.12557621
#> [2,] -0.11247063  1.00000000  0.07607558  0.19724379  0.08485137
#> [3,] -0.02124780  0.07607558  1.00000000 -0.10948813 -0.04034978
#> [4,]  0.06779336  0.19724379 -0.10948813  1.00000000  0.08710798
#> [5,] -0.12557621  0.08485137 -0.04034978  0.08710798  1.00000000
#> 

apply(x1, 2, function(x) moran.mc(x, listw = lw1, nsim = 2)$statistic)
#> [1]  0.08901257 -0.01819323  0.03040565 -0.02081023  0.01448972
apply(x1.6[[1]], 2, function(x) moran.mc(x, listw = lw1, nsim = 2)$statistic)
#> [1]  0.08901257 -0.01819323  0.03040565 -0.02081023  0.01448972