MCA {codep} | R Documentation |
Multiple-descriptors, Multiscale Codependence Analysis
Description
Class, Functions, and methods to perform Multiscale Codependence Analysis (MCA)
Usage
MCA(Y, X, emobj)
test.cdp(object, alpha = 0.05, max.step, response.tests = TRUE)
permute.cdp(object, permute, alpha = 0.05, max.step, response.tests = TRUE)
parPermute.cdp(
object,
permute,
alpha = 0.05,
max.step,
response.tests = TRUE,
nnode,
seeds,
verbose = TRUE,
...
)
Arguments
Y |
A numeric matrix or vector containing the response variable(s). |
X |
A numeric matrix or vector containing the explanatory variable(s). |
emobj |
A eigenmap-class object. |
object |
A cdp-class object. |
alpha |
The type I (alpha) error threshold used by the testing procedure. |
max.step |
The maximum number of steps to perform when testing for statistical significance. |
response.tests |
A boolean specifying whether to test individual response variables. |
permute |
The number of random permutations used for testing. When
omitted, the number of permutations is calculated using function
|
nnode |
The number of parallel computation nodes. |
seeds |
Seeds for computation nodes' random number generators when using parallel computation during the permutation test. |
verbose |
Whether to return user notifications. |
... |
Parameters to be passed to function |
Details
Multiscale Codependence Analysis (MCA) allows to calculate correlation-like (i.e.codependence) coefficients between two variables with respect to structuring variables (Moran's eigenvector maps). The purpose of this function is limited to parameter fitting.
Test procedures are handled through test.cdp
(parametric testing) or
permute.cdp
(permutation testing). Moreover, methods are provided for
printing (print.cdp
), displaying a summary of the tests
(summary.cdp
), plotting results (plot.cdp
), calculating
fitted (fitted.cdp
) and residuals values (redisuals.cdp
), and
making predictions (predict.cdp
).
It is noteworthy that the test procedure used by MCA
deviates from the
standard R workflow since intermediate testing functions (test.cdp
and
permute.cdp
) need first to be called before any testing be performed.
Function parPermute.cdp
allows the user to spread the number of
permutation on many computation nodes. It relies on package parallel.
Omitting argument nnode
lets function parallel::detectCores
specify the number of node. Similarly, omitting parameter seeds
lets
the function define the seeds as a set of values drawn from a uniform random
distribution between with minimum value -.Machine$integer.max
and
maximum value .Machine$integer.max
.
Value
A cdp-class object.
Functions
-
MCA()
: Main function to compute the multiscale codependence analysis -
test.cdp()
: Parametric statistical testing for multiscale codependence analysis -
permute.cdp()
: Permutation testing for multiscale codependence analysis. -
parPermute.cdp()
: Permutation testing for multiscale codependence analysis using parallel processing.
Author(s)
Guillaume Guenard and Pierre Legendre, Bertrand Pages Maintainer: Guillaume Guenard <guillaume.guenard@gmail.com>
References
Guénard, G., Legendre, P., Boisclair, D., and Bilodeau, M. 2010. Multiscale codependence analysis: an integrated approach to analyse relationships across scales. Ecology 91: 2952-2964
Guénard, G. Legendre, P. 2018. Bringing multivariate support to multiscale codependence analysis: Assessing the drivers of community structure across spatial scales. Meth. Ecol. Evol. 9: 292-304
Examples
### Example 1: St. Marguerite River Salmon Transect
data(salmon)
## Converting the data from data frames to to matrices:
Abundance <- log1p(as.matrix(salmon[,"Abundance",drop = FALSE]))
Environ <- as.matrix(salmon[,3L:5])
## Creating a spatial eigenvector map:
map1 <- eigenmap(x = salmon[,"Position"], weighting = wf.binary,
boundaries = c(0,20))
## Case of a single descriptor:
mca1 <- MCA(Y = Abundance, X = Environ[,"Substrate",drop = FALSE],
emobj = map1)
mca1
mca1_partest <- test.cdp(mca1)
mca1_partest
summary(mca1_partest)
par(mar = c(6,4,2,4))
plot(mca1_partest, las = 3, lwd=2)
mca1_pertest <- permute.cdp(mca1)
## Not run:
## or:
mca1_pertest <- parPermute.cdp(mca1, permute = 999999)
## End(Not run)
mca1_pertest
summary(mca1_pertest)
plot(mca1_pertest, las = 3)
mca1_pertest$UpYXcb$C # Array containing the codependence coefficients
## With all descriptors at once:
mca2 <- MCA(Y = log1p(as.matrix(salmon[,"Abundance",drop = FALSE])),
X = as.matrix(salmon[,3L:5]), emobj = map1)
mca2
mca2_partest <- test.cdp(mca2)
mca2_partest
summary(mca2_partest)
par(mar = c(6,4,2,4))
plot(mca2_partest, las = 3, lwd=2)
mca2_pertest <- permute.cdp(mca2)
## Not run:
## or:
mca2_pertest <- parPermute.cdp(mca2, permute = 999999)
## End(Not run)
mca2_pertest
summary(mca2_pertest)
plot(mca2_pertest, las = 3, lwd=2)
mca2_pertest$UpYXcb$C # Array containing the codependence coefficients
mca2_pertest$UpYXcb$C[,1L,] # now turned into a matrix.
### Example 2: Doubs Fish Community Transect
data(Doubs)
## Sites with no fish observed are excluded:
excl <- which(rowSums(Doubs.fish) == 0)
## Creating a spatial eigenvector map:
map2 <- eigenmap(x = Doubs.geo[-excl,"DFS"])
## The eigenvalues are in map2$lambda, the MEM eigenvectors in matrix map2$U
## MCA with multivariate response data analyzed on the basis of the Hellinger
## distance:
Y <- LGTransforms(Doubs.fish[-excl,],"Hellinger")
mca3 <- MCA(Y = Y, X=Doubs.env[-excl,], emobj = map2)
mca3_pertest <- permute.cdp(mca3)
## Not run:
## or:
mca3_pertest <- parPermute.cdp(mca3, permute = 999999)
## End(Not run)
mca3_pertest
summary(mca3_pertest)
par(mar = c(6,4,2,4))
plot(mca3_pertest, las = 2, lwd=2)
## Array containing all the codependence coefficients:
mca3_pertest$UpYXcb$C
## Display the results along the transect
spmeans <- colMeans(Y)
pca1 <- svd(Y - rep(spmeans, each=nrow(Y)))
par(mar = c(5,5,2,5) + 0.1)
plot(y = pca1$u[,1L], x = Doubs.geo[-excl,"DFS"], pch = 21L, bg = "red",
ylab = "PCA1 loadings", xlab = "Distance from river source (km)")
## A regular transect of sites from 0 to 450 (km) spaced by 1 km intervals
## (451 sites in total). It is used for plotting spatially-explicit
## predictions.
x <- seq(0,450,1)
newdists <- matrix(NA, length(x), nrow(Doubs.geo[-excl,]))
for(i in 1L:nrow(newdists))
newdists[i,] <- abs(Doubs.geo[-excl,"DFS"] - x[i])
## Calculating predictions for the regular transect under the same set of
## environmental conditions from which the codependence model was built.
prd1 <- predict(mca3_pertest,
newdata = list(target = eigenmap.score(map2, newdists)))
## Projection of the predicted species abundance on pca1:
Uprd1 <-
(prd1 - rep(spmeans, each = nrow(prd1))) %*%
pca1$v %*% diag(pca1$d^-1)
lines(y = Uprd1[,1L], x = x, col=2, lty = 1)
## Projection of the predicted species abundance on pca2:
plot(y = pca1$u[,2L], x = Doubs.geo[-excl,"DFS"], pch = 21L, bg = "red",
ylab = "PCA2 loadings", xlab = "Distance from river source (km)")
lines(y = Uprd1[,2L], x = x, col = 2, lty = 1)
## Displaying only the observed and predicted abundance for Brown Trout.
par(new = TRUE)
plot(y = Y[,"TRU"], Doubs.geo[-excl,"DFS"], pch = 21L,
bg = "green", ylab = "", xlab = "", new = FALSE, axes = FALSE)
axis(4)
lines(y = prd1[,"TRU"], x = x, col = 3)
mtext(side = 4, "sqrt(Brown trout rel. abundance)", line = 2.5)
### Example 3: Borcard et al. Oribatid Mite
## Testing the (2-dimensional) spatial codependence between the Oribatid Mite
## community structure and environmental variables, while displaying the
## total effects of the significant variables on the community structure
## (i.e., its first principal component).
data(mite)
map3 <- eigenmap(x = mite.geo)
Y <- LGTransforms(mite.species, "Hellinger")
## Organize the environmental variables
mca4 <- MCA(Y = Y, X = mite.env, emobj = map3)
mca4_partest <- test.cdp(mca4, response.tests = FALSE)
summary(mca4_partest)
plot(mca4_partest, las = 2, lwd = 2)
plot(mca4_partest, col = rainbow(1200)[1L:1000], las = 3, lwd = 4,
main = "Codependence diagram", col.signif = "white")
## Making a regular point grid for plotting the spatially-explicit
## predictions:
rng <- list(
x = seq(min(mite.geo[,"x"]) - 0.1, max(mite.geo[,"x"]) + 0.1, 0.05),
y = seq(min(mite.geo[,"y"]) - 0.1, max(mite.geo[,"y"]) + 0.1, 0.05))
grid <- cbind(x = rep(rng[["x"]], length(rng[["y"]])),
y = rep(rng[["y"]], each = length(rng[["x"]])))
newdists <- matrix(NA, nrow(grid), nrow(mite.geo))
for(i in 1L:nrow(grid)) {
newdists[i,] <- ((mite.geo[,"x"] - grid[i,"x"])^2 +
(mite.geo[,"y"] - grid[i,"y"])^2)^0.5
}
spmeans <- colMeans(Y)
pca2 <- svd(Y - rep(spmeans, each = nrow(Y)))
prd2 <- predict(mca4_partest,
newdata = list(target = eigenmap.score(map3, newdists)))
Uprd2 <-
(prd2 - rep(spmeans, each = nrow(prd2))) %*%
pca2$v %*% diag(pca2$d^-1)
## Printing the response variable (first principal component of the mite
## community structure).
prmat <- Uprd2[,1L]
dim(prmat) <- c(length(rng$x), length(rng$y))
zlim <- c(min(min(prmat), min(pca2$u[,1L])), max(max(prmat),
max(pca2$u[,1L])))
image(z = prmat, x = rng$x, y = rng$y, asp = 1, zlim = zlim,
col = rainbow(1200L)[1L:1000], ylab = "y", xlab = "x")
points(
x=mite.geo[,"x"], y=mite.geo[,"y"], pch=21L,
bg = rainbow(1200L)[round(1+(999*(pca2$u[,1L]-zlim[1L])/(zlim[2L]-zlim[1L])),0)])