Title: | Analysis of Habitat Selection by Animals |
---|---|
Description: | A collection of tools for the analysis of habitat selection. |
Authors: | Clement Calenge [aut, cre], Mathieu Basille [ctb] |
Maintainer: | Clement Calenge <[email protected]> |
License: | GPL (>= 2) |
Version: | 0.3.18 |
Built: | 2025-02-03 04:54:48 UTC |
Source: | https://github.com/clementcalenge/adehabitaths |
This data set contains the relocations of 198 chamois groups in the Bauges mountains, as well as maps of 7 environmental variables on the study area.
data(bauges)
data(bauges)
This dataset contains a subsample of the data collected by volunteers and professionals working in various French wildlife and forest management, from 1994 to 2004 in the wildlife reserve of "les Bauges" (French Alps). Note that both the maps and the relocations have been slightly destroyed to preserve copyright on the data.
Daniel Maillard, Office National de la chasse et de la faune sauvage, 95 rue Pierre Flourens, 34000 Montpellier, France
data(bauges) mimage(bauges$map) image(bauges$map, 1) points(bauges$locs)
data(bauges) mimage(bauges$map) image(bauges$map, 1) points(bauges$locs)
This data set describes the habitat use and availability for 6 bighorn sheeps monitored by radio-tracking (Arnett et al. 1989, in Manly et al., 2003, p. 67-74). 10 habitat types are considered.
data(bighorn)
data(bighorn)
The object bighorn
is a list, with the following components:
used
the number of resource units used by each animal (in rows) in each habitat category (in columns).
availTrue
the availability of each habitat category.
availEstimated
a sample of available resource units in each habitat category.
Manly, B.F.J., McDonald, L.L., Thomas, D.L., McDonald, T.L. & Erickson, W.P. (2003) Resource selection by animals - Statistical design and Analysis for field studies. Second edition. London: Kluwer academic publishers.
biv.plot
displays a bivariate plot.
biv.test
displays the results of a bivariate randomisation test.
biv.plot(dfxy, br = 10, points = TRUE, density = TRUE, kernel = TRUE, o.include = FALSE, pch, cex, col, h, sub, side = c("top", "bottom", "none"), ...) biv.test(dfxy, point, br = 10, points = TRUE, density = TRUE, kernel = TRUE, o.include = FALSE, pch, cex, col, Pcol, h, sub, side = c("top", "bottom", "none"), ...)
biv.plot(dfxy, br = 10, points = TRUE, density = TRUE, kernel = TRUE, o.include = FALSE, pch, cex, col, h, sub, side = c("top", "bottom", "none"), ...) biv.test(dfxy, point, br = 10, points = TRUE, density = TRUE, kernel = TRUE, o.include = FALSE, pch, cex, col, Pcol, h, sub, side = c("top", "bottom", "none"), ...)
dfxy |
a data frame with N lines (couples of values) and two columns |
br |
a parameter used to define the numbers of breaks of the histograms. A larger value leads to a larger number of breaks |
points |
logical. Whether the points should be displayed |
density |
logical. Whether the kernel density estimation should be displayed for the marginal histograms |
kernel |
logical. Whether the kernel density estimation should be displayed for the bivariate plot |
o.include |
logical. If |
pch |
plotting "character", i.e., symbol to use for the points. (see
|
cex |
character expansion for the points |
col |
color code or name for the points, see |
h |
vector of bandwidths for x and y directions, used in the
function |
sub |
a character string to be inserted in the plot as a title |
side |
if |
point |
a vector of length 2, representing the observation to be compared with the simulated values of the randomisation test |
Pcol |
color code or name for the observed point |
... |
further arguments passed to or from other methods |
biv.test
is used to display the results of a bivariate
randomisation test. An example of use of the function is provided in
the function niche.test
.
The x-axis of the main window corresponds to the first column of
dfxy
; the y-axis corresponds to the second column. Kernel
density is estimated to indicate the contours of the distribution of
randomised values. The two marginal histograms correspond to the
univariate tests on each axis, for which the p-values are computed with
as.randtest
(package ade4, one-sided tests).
biv.plot
and biv.test
use the function kde2d
of the
package MASS
.
Mathieu Basille [email protected]
as.randtest
(package ade4)
x = rnorm(1000,2) y = 2*x+rnorm(1000,2) dfxy = data.frame(x, y) biv.plot(dfxy) biv.plot(dfxy, points=FALSE, col="lightblue", br=20) p = c(3, 4) biv.test(dfxy, p) biv.test(dfxy, p, points=FALSE, Pcol="darkred", col="lightblue", br=20)
x = rnorm(1000,2) y = 2*x+rnorm(1000,2) dfxy = data.frame(x, y) biv.plot(dfxy) biv.plot(dfxy, points=FALSE, col="lightblue", br=20) p = c(3, 4) biv.test(dfxy, p) biv.test(dfxy, p, points=FALSE, Pcol="darkred", col="lightblue", br=20)
This function performs a canonical OMI analysis (outlying mean index).
canomi(dudiX, Y, scannf = TRUE, nf = 2) ## S3 method for class 'canomi' print(x, ...) ## S3 method for class 'canomi' plot(x, xax = 1, yax = 2, ...)
canomi(dudiX, Y, scannf = TRUE, nf = 2) ## S3 method for class 'canomi' print(x, ...) ## S3 method for class 'canomi' plot(x, xax = 1, yax = 2, ...)
dudiX |
an object of class |
Y |
a a data frame Resource units-animals according to
|
scannf |
a logical value indicating whether the eigenvalues bar plot should be displayed |
nf |
if scannf FALSE, an integer indicating the number of kept axes |
x |
an object of class |
xax |
the number of the x-axis |
yax |
the number of the x-axis |
... |
further arguments passed to or from other methods |
The canonical OMI analysis is similar to the function niche
,
from the package ade4
. The principle of this analysis is the
following. A set of N resource units (RUs) are available to the
animals of the study. Each resource unit is described by P
environmental variables. Therefore, the N resource units define a
cloud of N points in the P-dimensionnal space defined by the P
variables. We call this space "ecological space".
Moreover, the use of the N resource units is known (or sampled) for
a sample of K animals (e.g., using radio-tracking). These utilization
weights for each RU (rows) and each animal (column) define a table
Y
. For a given animal, the set of resource units used define
the "niche" of the animal. The vector connecting the centroid (mean)
of the available RUs to the centroid of the RUs used by this animal is
named "marginality vector" (and its squared length is named
"marginality" or "outlying mean index").
The canomi
first distorts the ecological space so that the
available resource units take a standard spherical shape (by first
performing a principal component analysis). Then, in this distorted
space, a non-centred principal component analysis of the marginality
vectors is performed. The canonical OMI analysis finds the directions
in the distorted ecological space where the marginality is, in
average, the largest.
canomi
returns a list of the class canomi
, with the
following components:
call |
original call. |
rank |
an integer indicating the rank of the studied matrix |
nf |
an integer indicating the number of kept axes |
eig |
a vector with all the eigenvalues of the analysis. |
tab |
a data frame with n rows (n animals) and p columns (p environmental variables). |
li |
animals coordinates, data frame with n rows and nf columns. |
l1 |
animals normed coordinates, data frame with n rows and nf columns. |
c1 |
column scores, data frame with p rows and nf columns. |
cor |
the correlation between the canomi axes and the original variables |
ls |
a data frame with the resource units coordinates |
cm |
The variables metric used in the analysis (e.g.
|
as |
a data frame with the axis upon niche axis |
Clement Calenge [email protected]
Chessel, D. 2006. Calcul de l'outlier mean index. Consultation statistique avec le logiciel R.
dudi
for class dudi
,
niche
for classical OMI analysis
## The data data(puech) locs <- puech$relocations maps <- puech$maps ## the maps mimage(maps) ## the relocations of the wild boar: image(maps) points(locs, col=as.numeric(slot(locs, "data")[,1])) ## count the number of relocations ## in each pixel of the maps cp <- count.points(locs, maps) ## gets the data: dfavail <- slot(maps, "data") dfused <- slot(cp, "data") ## a preliminary principal component analysis of the data: dud <- dudi.pca(dfavail, scannf=FALSE) ## The analysis: nic <- canomi(dud, dfused, scannf=FALSE) nic ## Plot the results: plot(nic)
## The data data(puech) locs <- puech$relocations maps <- puech$maps ## the maps mimage(maps) ## the relocations of the wild boar: image(maps) points(locs, col=as.numeric(slot(locs, "data")[,1])) ## count the number of relocations ## in each pixel of the maps cp <- count.points(locs, maps) ## gets the data: dfavail <- slot(maps, "data") dfused <- slot(cp, "data") ## a preliminary principal component analysis of the data: dud <- dudi.pca(dfavail, scannf=FALSE) ## The analysis: nic <- canomi(dud, dfused, scannf=FALSE) nic ## Plot the results: plot(nic)
This data set describes the habitat use and availability by the chamois of the Chartreuse mountains (Isere, France), in 1992 and 1997. These data have been gathered during the hunting season (Fall).
data(chamois)
data(chamois)
The object chamois
is a list containing the following
components:
locs
is a data frame containing the x and y coordinates
of 198 chamois groups.
map
is a map of class kasc
describing the
vegetation (Forest or Open areas), the distance from the ecotone
Open/Forest and the slopes on the area.
Federation Departementale des Chasseurs de l'Isere, 65 av Jean Jaures, 38320 Eybens. France.
compana
performs a classical compositional analysis of habitat
use (Aebischer et al., 1993).
compana(used, avail, test = c("randomisation", "parametric"), rnv = 0.01, nrep = 500, alpha = 0.1)
compana(used, avail, test = c("randomisation", "parametric"), rnv = 0.01, nrep = 500, alpha = 0.1)
used |
a matrix or a data frame describing the percentage of use of habitats (in columns) by animals (in rows). |
avail |
a matrix or a data frame describing the percentage of availability of habitats (in columns) by animals (in rows). |
test |
a character string. If |
rnv |
the number replacing the 0 values occurring in the
matrix |
nrep |
the number of repetitions in the randomisation tests. |
alpha |
the alpha level for the tests. |
The compositional analysis of habitat use has been recommended by Aebischer et al. (1993) for the analysis of habitat selection by several animals, when the resources are defined by several categories (e.g. vegetation types).
This analysis is carried out in two steps: first the significance of habitat selection is tested (using a Wilks lambda). Then, a ranking matrix is built, indicating whether the habitat type in row is significantly used more or less than the habitat type in column. When this analysis is performed on radio-tracking data, Aebischer et al. recommend to study habitat selection at two levels: (i) selection of the home range within the study area, and (ii) selection of the relocations within the home range. The first level is termed second-order habitat selection on Johnson's scale (1980), and the second one, third-order habitat selection.
When zero values are found in the matrix of used habitats, they are replaced by a small value (by default, 0.01), according to the recommendations of Aebischer et al. (1993).
When zero values are found in the matrix of available habitats, the
function compana
uses the procedure termed
"weighted mean lambda" described in Aebischer et al.
(1993: Appendix 2), instead of the usual lambda (see examples). Zero
values can be found in the matrix of available habitats when the
third-order habitat selection is under focus. In this case, it may
occur that some habitat types are available to some animals and not to
the others.
Note that this method rely on the following hypotheses: (i)
independence between animals, and (ii) all animals are selecting
habitat in the same way (in addition to "traditional" hypotheses in
these kinds of studies: no territoriality, all animals having equal
access to all available resource units, etc.). The function
eisera
can be used as a preliminary to identify whether this is
indeed the case (see examples).
Returns a list of the class compana
:
used |
the matrix of used habitats |
avail |
the matrix of available habitats |
type.test |
a character string. Either |
test |
the results of the test of habitat selection |
rm |
the ranking matrix: a square matrix with |
rmnb |
the matrix containing the number of animals used to
perform the tests in |
rank |
the rank of the habitat types. It is equal to the number
of |
rmv |
the matrix of statistics used to build |
profile |
the profile of preferences: resource types are sorted so that the left type is the most preferred and the right type is the most avoided. Habitats for which the intensity of habitat selection is similar (no significant difference) are connected by a line. |
In the examples below, the results differ from those published in Aebischer et al. (squirrel example, selection of the relocations within the home range). In fact, there has been a confusion in the column names in the paper. Actually, Aebischer (pers. com.) indicated that the ranking matrix given in this example is correct.
Clement Calenge [email protected]
Aebischer, N. J. and Robertson, P. A. (1992) Practical aspects of compositional analysis as applied to pheasant habitat utilisation. pp. 285–293 In: Priede, G. and Swift, S. M. Wildlife telemetry, remote monitoring and tracking of animals.
Aebischer, N. J., Robertson, P. A. and Kenward, R. E. (1993) Compositional analysis of habitat use from animal radiotracking data. Ecology, 74, 1313–1325.
Johnson, D. H. (1980) The comparison of usage and availability measurements for evaluating resource preference. Ecology, 61, 65–71.
eisera
to perform an eigenanalysis of selection
ratios, preliminary to the use of compositional analysis.
## The examples presented here ## are the same as those presented in ## the paper of Aebischer et al. (1993) ############################# ## Pheasant dataset: first ## example in Aebischer et al. data(pheasant) ## Second order habitat selection ## Selection of home range within the ## study area (example of parametric test) pheana2 <- compana(pheasant$mcp, pheasant$studyarea, test = "parametric") pheana2 ## The ranking matrix: print(pheana2$rm, quote = FALSE) ## Third order habitat selection ## (relocation within home range) ## We remove the first pheasant of the analysis ## (as in the paper of Aebischer et al.) ## before the analysis pheana3 <- compana(pheasant$locs[-1,], pheasant$mcp[-1,c(1,2,4)]) pheana3 ## The ranking matrix: print(pheana3$rm, quote = FALSE) ############################# ## Squirrel data set: second ## example in Aebischer et al. data(squirrel) ## Second order habitat selection ## Selection of home range within the ## study area squiana2 <- compana(squirrel$mcp, squirrel$studyarea) squiana2 ## The ranking matrix: print(squiana2$rm, quote = FALSE) ## However, note that here, the hypothesis of identical use ## on which this analysis relies is likely to be false. ## Indeed, an eisera indicates: us <- round(30 * squirrel$locs / 100) av <- squirrel$studyarea ii <- eisera(us, av, scannf = FALSE) scatter(ii, grid = FALSE, clab = 0.7) ## There are clearly two groups of animals. In such cases, ## compositional analysis is to be avoided in this case. ## Third order habitat selection ## (relocation within home range) ## We remove the second column ## (as in the paper of Aebischer et al.) squiana3 <- compana(squirrel$locs[,-2], squirrel$mcp[,-2]) squiana3 ## The ranking matrix: print(squiana3$rm, quote = FALSE)
## The examples presented here ## are the same as those presented in ## the paper of Aebischer et al. (1993) ############################# ## Pheasant dataset: first ## example in Aebischer et al. data(pheasant) ## Second order habitat selection ## Selection of home range within the ## study area (example of parametric test) pheana2 <- compana(pheasant$mcp, pheasant$studyarea, test = "parametric") pheana2 ## The ranking matrix: print(pheana2$rm, quote = FALSE) ## Third order habitat selection ## (relocation within home range) ## We remove the first pheasant of the analysis ## (as in the paper of Aebischer et al.) ## before the analysis pheana3 <- compana(pheasant$locs[-1,], pheasant$mcp[-1,c(1,2,4)]) pheana3 ## The ranking matrix: print(pheana3$rm, quote = FALSE) ############################# ## Squirrel data set: second ## example in Aebischer et al. data(squirrel) ## Second order habitat selection ## Selection of home range within the ## study area squiana2 <- compana(squirrel$mcp, squirrel$studyarea) squiana2 ## The ranking matrix: print(squiana2$rm, quote = FALSE) ## However, note that here, the hypothesis of identical use ## on which this analysis relies is likely to be false. ## Indeed, an eisera indicates: us <- round(30 * squirrel$locs / 100) av <- squirrel$studyarea ii <- eisera(us, av, scannf = FALSE) scatter(ii, grid = FALSE, clab = 0.7) ## There are clearly two groups of animals. In such cases, ## compositional analysis is to be avoided in this case. ## Third order habitat selection ## (relocation within home range) ## We remove the second column ## (as in the paper of Aebischer et al.) squiana3 <- compana(squirrel$locs[,-2], squirrel$mcp[,-2]) squiana3 ## The ranking matrix: print(squiana3$rm, quote = FALSE)
domain
uses the DOMAIN algorithm to estimate the potential
distribution of a species based on a list of species occurrences and on
maps of the area.
domain(x, pts, type = c("value", "potential"), thresh = 0.95)
domain(x, pts, type = c("value", "potential"), thresh = 0.95)
x |
an object of class |
pts |
a data frame giving the x and y coordinates of the species occurrences. |
type |
a character string. The |
thresh |
if |
This function implements the DOMAIN algorithm described in Carpenter et al. (1993).
Returns an object of class SpatialPixelsDataFrame
.
domain
is restricted to maps
containing only numerical variables (i.e. no factors).
Clement Calenge [email protected]
Carpenter, G., Gillison, A.N. and Winter, J. (1993) DOMAIN: a flexible modelling procedure for mapping potential distributions of plants and animals. Biodiversity and conservation, 2, 667–680.
## Preparation of the data data(lynxjura) map <- lynxjura$map pts <- lynxjura$locs ## View of the data image(map) title(main="Elevation") points(pts, pch = 3) ## Estimation of habitat suitability map hsm <- domain(map, pts) image(hsm, col = grey((1:256)/256)) contour(hsm, add = TRUE) ## Lighter areas are the most used areas ## Potential distribution hsm <- domain(map, pts, type = "potential", thresh = 0.98) image(hsm, col = "orange") title(main = "Habitat suitability map") points(pts, pch = 3)
## Preparation of the data data(lynxjura) map <- lynxjura$map pts <- lynxjura$locs ## View of the data image(map) title(main="Elevation") points(pts, pch = 3) ## Estimation of habitat suitability map hsm <- domain(map, pts) image(hsm, col = grey((1:256)/256)) contour(hsm, add = TRUE) ## Lighter areas are the most used areas ## Potential distribution hsm <- domain(map, pts, type = "potential", thresh = 0.98) image(hsm, col = "orange") title(main = "Habitat suitability map") points(pts, pch = 3)
dunnfa
performs a factorial decomposition of the Mahalanobis
distances in habitat selection studies (see details).
dunnfa(dudi, pr, scannf = TRUE, nf = 2) ## S3 method for class 'dunnfa' print(x, ...)
dunnfa(dudi, pr, scannf = TRUE, nf = 2) ## S3 method for class 'dunnfa' print(x, ...)
dudi |
an object of class |
pr |
a vector giving the utilization weights associated to each unit |
scannf |
logical. Whether the eigenvalues barplot should be displayed |
nf |
an integer indicating the number of kept axes |
x |
an object of class |
... |
additional arguments to be passed to the function
|
This analysis is in essence very similar to the MADIFA (see
?madifa
). The Mahalanobis distances are often used in the
context of niche-environment studies (Clark et al. 1993, see the
function mahasuhab
). Each resource unit takes a value on a
set of environmental variables. Each environmental variable defines a
dimension in a multidimensionnal space, namely the ecological space.
A set of points (resource units) describes what is available to the
species. For each point, a "utilization weight" measures the
intensity of use of the point by the species. The set of points for
which the utilization weight is greater than zero defines the
"niche". The Mahalanobis distance between any resource unit in this
space (e.g. the point defined by the values of environmental variables
in a pixel of a raster map) and the centroid of the niche (the
distribution of used resource units) can be used to give a value of
eccentricity to this point.
For a given distribution of available resource units, for which a
measure of Mahalanobis distances is desired, the MADIFA (MAhalanobis
DIstances Factor Analysis) partitions the ecological space into a set
of axes, so that the first axes maximises the average proportion of
their squared Mahalanobis distances. James Dunn (formerly University
of Arkansas) proposed the analysis programmed in the function
dunnfa
, as an alternative to the MADIFA (unpublished
results). This analysis is closely related to both the ENFA
(Ecological niche factor analysis, Hirzel et al. 2002) and the
MADIFA.
The analysis proposed by James Dunn searches, in the multidimensional space defined by environmental variables, synthesis variables which maximise the ratio (variance of the scores of available resource units) / (variance of the scores of used resource units). This ratio is sometimes called "specialization" in the ecological literature (Hirzel et al. 2002). It is therefore very similar to the ENFA (which also maximises the specialization), except that the factorial axes returned by this analysis are not required to be *orthogonal to the marginality axis*.
James Dunn demonstrated that this analysis also partitions the Mahalanobis distances into uncorrelated axes, which makes it similar to the MADIFA (the difference is that the MADIFA maximises the mean squared Mahalanobis distances on the first axes, whereas the DUNNFA maximises the specialization on the first axes). Therefore, as for the MADIFA, the DUNNFA can be used to build reduced rank habitat suitability map.
Note that although this analysis could theoretically be used with all kinds of variables, it it currently implemented only for numeric variables.
dunnfa
returns a list of class dunnfa
containing the
following components:
call |
original call. |
tab |
a data frame with n rows and p columns (original data frame centered by column for the uniform weighting). |
pr |
a vector of length n containing the number of points in each pixel of the map. |
nf |
the number of kept axes. |
eig |
a vector with all the eigenvalues of the analysis. |
liA |
row coordinates (centering on the centroid of the cloud of available points), data frame with n rows and nf columns. |
liU |
row coordinates (centering on the centroid of the cloud of available points), data frame with p rows and nf columns. |
mahasu |
a vector of length n containing the reduced-rank squared Mahalanobis distances for the n units. |
co |
column (environmental variables) coordinates, a data frame with p rows and nf columns |
cor |
the correlation between the DUNNFA axes and the original variable |
This analysis was developed by James Dunn during an e-mail discussion on the MADIFA, and is still unpublished work. Implemented in adehabitatHS with his autorization.
Clement Calenge [email protected]
Clark, J.D., Dunn, J.E. and Smith, K.G. (1993) A multivariate model of
female black bear habitat use for a geographic information
system. Journal of Wildlife Management, 57, 519–526.
Hirzel, A.H., Hausser, J., Chessel, D. & Perrin, N. (2002) Ecological-niche factor analysis: How to compute habitat-suitability maps without absence data? Ecology, 83, 2027–2036.
Calenge, C., Darmon, G., Basille, M., Loison, A. and Jullien
J.M. (2008) The factorial decomposition of the Mahalanobis distances
in habitat selection studies. Ecology, 89, 555–566.
madifa
, enfa
and
gnesfa
for related methods. mahasuhab
for
details about the Mahalanobis distances.
## Not run: data(bauges) map <- bauges$map locs <- bauges$loc ## We prepare the data for the analysis tab <- slot(map, "data") pr <- slot(count.points(locs, map), "data")[,1] ## We then perform the PCA before the analysis pc <- dudi.pca(tab, scannf = FALSE) (dun <- dunnfa(pc, pr, nf=2, scannf = FALSE)) ## We should keep one axis: barplot(dun$eig) ## The correlation of the variables with the first two axes: s.arrow(dun$cor) ## A factorial map of the niche (centering on the available points) scatterniche(dun$liA, dun$pr, pts=TRUE) ## a map of the reduced rank Mahalanobis distances ## (here, with one axis) dun2 <- dunnfa(pc, pr, nf=1, scannf = FALSE) df <- data.frame(MD=dun2$mahasu) coordinates(df) <- coordinates(map) gridded(df) <- TRUE image(df) ## Compute the specialization on the row scores of ## the analysis: apply(dun$liA, 2, function(x) { varav <- sum((x - mean(x))^2) / length(x) meanus <- sum(dun$pr*x)/sum(dun$pr) varus <- sum(dun$pr * (x - meanus )^2)/sum(dun$pr) return(varav/varus) }) ## The eigenvalues: dun$eig ## End(Not run)
## Not run: data(bauges) map <- bauges$map locs <- bauges$loc ## We prepare the data for the analysis tab <- slot(map, "data") pr <- slot(count.points(locs, map), "data")[,1] ## We then perform the PCA before the analysis pc <- dudi.pca(tab, scannf = FALSE) (dun <- dunnfa(pc, pr, nf=2, scannf = FALSE)) ## We should keep one axis: barplot(dun$eig) ## The correlation of the variables with the first two axes: s.arrow(dun$cor) ## A factorial map of the niche (centering on the available points) scatterniche(dun$liA, dun$pr, pts=TRUE) ## a map of the reduced rank Mahalanobis distances ## (here, with one axis) dun2 <- dunnfa(pc, pr, nf=1, scannf = FALSE) df <- data.frame(MD=dun2$mahasu) coordinates(df) <- coordinates(map) gridded(df) <- TRUE image(df) ## Compute the specialization on the row scores of ## the analysis: apply(dun$liA, 2, function(x) { varav <- sum((x - mean(x))^2) / length(x) meanus <- sum(dun$pr*x)/sum(dun$pr) varus <- sum(dun$pr * (x - meanus )^2)/sum(dun$pr) return(varav/varus) }) ## The eigenvalues: dun$eig ## End(Not run)
Performs an eigenanalysis of selection ratios.
eisera(used, available, scannf = TRUE, nf = 2) ## S3 method for class 'esr' print(x, ...) ## S3 method for class 'esr' scatter(x, xax = 1, yax = 2, csub = 1, possub = "bottomleft", ...)
eisera(used, available, scannf = TRUE, nf = 2) ## S3 method for class 'esr' print(x, ...) ## S3 method for class 'esr' scatter(x, xax = 1, yax = 2, csub = 1, possub = "bottomleft", ...)
used |
a data frame containing the *number* of relocations of each animal (rows) in each habitat type (columns) |
available |
a data frame containing the *proportion* of availability of each habitat type (columns) to each animal (rows) |
scannf |
logical. Whether the eigenvalues bar plot should be displayed |
nf |
if |
x |
an object of class |
xax |
the column number for the x-axis |
yax |
the column number for the y-axis |
csub |
a character size for the legend, used with
|
possub |
a string of characters indicating the sub-title position ("topleft", "topright", "bottomleft", "bottomright") |
... |
further arguments passed to or from other methods |
The eigenanalysis of selection ratios has been developped to explore
habitat selection by animals monitored using radio-tracking, when
habitat is defined by several categories (e.g. several vegetation
types, see Calenge and Dufour 2006).
This analysis can be used for both designs II (same availability for all animals, e.g. selection of the home range within the study area) and designs III (different availability, e.g. selection of the sites within the home range). In the latter case, when some available proportions are equal to zero, the selection ratios are replaced by their expectation under random habitat use, following the recommendations of Calenge and Dufour (2006).
A list of class esr
and dudi
containing also:
available |
available proportions |
used |
number of relocations |
wij |
selection ratios |
Clement Calenge [email protected]
Calenge, C. and Dufour, A.B. (2006) Eigenanalysis of selection ratios from animal radio-tracking data. Ecology. 87, 2349–2355.
wi
for further information about the
selection ratios, compana
for compositional analysis.
########################################################### ########################################################### ### ### Example given in Calenge and Dufour 2006 (design II) data(squirrel) ## computation of the number of relocations in each habitat type ## from the data given by Aebischer et al. (1993). ## squirrel$locs give the percentage of relocations in each habitat ## type, and Aebischer et al. (1993) indicate that there are 30 ## relocations per animal. ## We therefore compute the number of relocations in each habitat type ## using: us <- round(30 * squirrel$locs / 100) ## Habitat availability av <- squirrel$studyarea ## Eigenanalysis of selection ratios ii <- eisera(us, av, scannf = FALSE) scatter(ii, grid = FALSE, clab = 0.7) ## The following graph may help the interpretation ## (see Calenge and Dufour 2006) data(squirreloc) locs <- squirreloc$locs are <- squirreloc$map ind <- levels(slot(locs, "data")$id) opar <- par(mfrow=n2mfrow(length(ind)), mar=c(0,0,2,0)) tmp <- lapply(1:length(ind), function(i) { plot(are, col = as.data.frame(are)[,2]) title(main = ind[i]) points(locs[slot(locs, "data")[,1]==ind[i],], pch=16, cex=1.5) box() }) plot(0,0, axes=FALSE, ty="n", xlim=c(-1,1), asp=1) legend(-0.8,0.8, unique(slot(are,"data")[,1]), fill=unique(slot(are,"data")[,2])) par(opar) ########################################################### ########################################################### ### ### Example of design III iii <- eisera(us, squirrel$mcp, scannf = FALSE) scatter(iii, grid = FALSE, clab = 0.7)
########################################################### ########################################################### ### ### Example given in Calenge and Dufour 2006 (design II) data(squirrel) ## computation of the number of relocations in each habitat type ## from the data given by Aebischer et al. (1993). ## squirrel$locs give the percentage of relocations in each habitat ## type, and Aebischer et al. (1993) indicate that there are 30 ## relocations per animal. ## We therefore compute the number of relocations in each habitat type ## using: us <- round(30 * squirrel$locs / 100) ## Habitat availability av <- squirrel$studyarea ## Eigenanalysis of selection ratios ii <- eisera(us, av, scannf = FALSE) scatter(ii, grid = FALSE, clab = 0.7) ## The following graph may help the interpretation ## (see Calenge and Dufour 2006) data(squirreloc) locs <- squirreloc$locs are <- squirreloc$map ind <- levels(slot(locs, "data")$id) opar <- par(mfrow=n2mfrow(length(ind)), mar=c(0,0,2,0)) tmp <- lapply(1:length(ind), function(i) { plot(are, col = as.data.frame(are)[,2]) title(main = ind[i]) points(locs[slot(locs, "data")[,1]==ind[i],], pch=16, cex=1.5) box() }) plot(0,0, axes=FALSE, ty="n", xlim=c(-1,1), asp=1) legend(-0.8,0.8, unique(slot(are,"data")[,1]), fill=unique(slot(are,"data")[,2])) par(opar) ########################################################### ########################################################### ### ### Example of design III iii <- eisera(us, squirrel$mcp, scannf = FALSE) scatter(iii, grid = FALSE, clab = 0.7)
enfa
performs an Ecological Niche Factor Analysis.
hist.enfa
draws histograms of the row scores or of the initial
variables of the ENFA.
enfa(dudi, pr, scannf = TRUE, nf = 1) ## S3 method for class 'enfa' hist(x, scores = TRUE, type = c("h", "l"), adjust = 1, Acol, Ucol, Aborder, Uborder, Alwd = 1, Ulwd = 1, ...)
enfa(dudi, pr, scannf = TRUE, nf = 1) ## S3 method for class 'enfa' hist(x, scores = TRUE, type = c("h", "l"), adjust = 1, Acol, Ucol, Aborder, Uborder, Alwd = 1, Ulwd = 1, ...)
dudi |
a duality diagram, object of class |
pr |
a vector giving the utilization weights associated to each unit |
scannf |
logical. Whether the eigenvalues barplot should be displayed |
nf |
an integer indicating the number of kept specialization axes |
x |
an object of class |
scores |
logical. If |
type |
what type of plot should be drawn. Possible types are: |
adjust |
if |
Acol |
if |
Ucol |
if |
Aborder |
color for the border of the histograms of the available pixels |
Uborder |
color for the border of the histograms of the used pixels |
Alwd |
if |
Ulwd |
if |
... |
further arguments passed to or from other methods |
The niche concept, as defined by Hutchinson (1957), considers the ecological niche of a species as an hypervolume in the multidimensional space defined by environmental variables, within which the populations of a species can persist. The Ecological Niche Factor Analysis (ENFA) has been developped by Hirzel et al. (2002) to analyse the position of the niche in the ecological space. Nicolas Perrin (1984) described the position of the niche in the n-dimensional space using two measures: the M-specialization (hereafter termed marginality) and the S-specialization (hereafter termed specialization). The marginality represents the squared distance of the niche barycentre from the mean available habitat. A large specialization corresponds to a narrow niche relative to the habitat conditions available to the species.
The ENFA first extracts an axis of marginality (vector from the average of available habitat conditions to the average used habitat conditions). Then the analysis extracts successives orthogonal axes (i.e. uncorrelated), which maximises the specialization of the species. The calculations used in the function are described in Hirzel et al. (2002).
The function enfa
can be used on both quantitative variables and
qualitative variables (though the interpretation of the results of the
ENFA for qualitative variables is still under research), provided that
the table containing the values of habitat variables (columns) for
each resource unit (rows) is correctly transformed
(e.g. column-centered and standardised for tables containing only
quantitative variables), and that appropriate column weights are given
(e.g. the sum of the weights for
the levels of a factor should be the same as the weight of one
quantitative variable). Therefore, the function enfa
requires
that a preliminary multivariate analysis is performed on the table
(using analysis of the family of duality diagram, e.g. principal
component analysis or Hill and Smith analysis). The object returned
by this preliminary analysis contains the appropriate weights and
transformation of the original data frame. For example,
the function dudi.mix
can be used first on the data.frame
containing the value of both quantitative (e.g. slope, elevation) and
qualitative habitat variables (e.g. vegetation) for each pixel of a
raster map. The result of this analysis can then be passed as argument
to the function enfa
(see examples below).
enfa
returns a list of class enfa
containing the
following components:
call |
original call. |
tab |
a data frame with n rows and p columns. |
pr |
a vector of length n containing the number of points in each pixel of the map. |
nf |
the number of kept specialization axes. |
m |
the marginality (squared length of the marginality vector). |
s |
a vector with all the eigenvalues of the analysis. |
lw |
row weights, a vector with n components. |
li |
row coordinates, data frame with n rows and nf columns. |
cw |
column weights, a vector with p components. |
co |
column coordinates, data frame with p rows and nf columns. |
mar |
coordinates of the marginality vector. |
Mathieu Basille [email protected]
Hutchinson, G.E. (1957) Concluding Remarks. Cold Spring Harbor Symposium on Quantitative Biology, 22: 415–427.
Perrin, N. (1984) Contribution a l'ecologie du genre Cepaea (Gastropoda) : Approche descriptive et experimentale de l'habitat et de la niche ecologique. These de Doctorat. Universite de Lausanne, Lausanne.
Hirzel, A.H., Hausser, J., Chessel, D. and Perrin, N. (2002) Ecological-niche factor analysis: How to compute habitat-suitability maps without absence data? Ecology, 83, 2027–2036.
Basille, M., Calenge, C., Marboutin, E., Andersen, R. and Gaillard, J.M. (2008) Assessing habitat selection using multivariate statistics: Some refinements of the ecological-niche factor analysis. Ecological Modelling, 211, 233–240.
niche
, kselect
for other types of
analysis of the niche, when several species are under studies,
and scatter.enfa
to have a
graphical display of objects of class
enfa
. See madifa
for another factorial analysis of
the ecological niche.
data(lynxjura) map <- lynxjura$map ## We keep only "wild" indices. locs <- lynxjura$locs locs <- locs[slot(locs, "data")[,2]!="D",] hist(map, type = "l") ## The variable artif is far from symetric ## We perform a square root transformation ## of this variable ## We therefore normalize the variable 'artif' slot(map,"data")[,4] <- sqrt(slot(map,"data")[,4]) hist(map, type = "l") ## We prepare the data for the ENFA tab <- slot(map, "data") pr <- slot(count.points(locs, map), "data")[,1] ## We then perform the PCA before the ENFA pc <- dudi.pca(tab, scannf = FALSE) ## The object 'pc' contains the transformed table (i.e. ## centered so that all columns have a mean of 0 ## and scaled so that all columns have a variance of 1 ## 'pc' also contains the weights of the habitat variables, ## and the weights of the pixels in the analysis (enfa1 <- enfa(pc, pr, scannf = FALSE)) hist(enfa1) hist(enfa1, scores = FALSE, type = "l") ## scatterplot scatter(enfa1) ## randomization test ## Not run: (renfa <- randtest(enfa1)) plot(renfa) ## End(Not run)
data(lynxjura) map <- lynxjura$map ## We keep only "wild" indices. locs <- lynxjura$locs locs <- locs[slot(locs, "data")[,2]!="D",] hist(map, type = "l") ## The variable artif is far from symetric ## We perform a square root transformation ## of this variable ## We therefore normalize the variable 'artif' slot(map,"data")[,4] <- sqrt(slot(map,"data")[,4]) hist(map, type = "l") ## We prepare the data for the ENFA tab <- slot(map, "data") pr <- slot(count.points(locs, map), "data")[,1] ## We then perform the PCA before the ENFA pc <- dudi.pca(tab, scannf = FALSE) ## The object 'pc' contains the transformed table (i.e. ## centered so that all columns have a mean of 0 ## and scaled so that all columns have a variance of 1 ## 'pc' also contains the weights of the habitat variables, ## and the weights of the pixels in the analysis (enfa1 <- enfa(pc, pr, scannf = FALSE)) hist(enfa1) hist(enfa1, scores = FALSE, type = "l") ## scatterplot scatter(enfa1) ## randomization test ## Not run: (renfa <- randtest(enfa1)) plot(renfa) ## End(Not run)
These functions implements the method described by Engen et al. to measure the preference of animals for habitat variables in habitat selection studies.
engen2008II(us, av, id, nsim = 500, nsimra = 500) engen2008I(us, av, nsimra=500) ## S3 method for class 'engenetalI' print(x, ...) ## S3 method for class 'engenetalII' print(x, ...)
engen2008II(us, av, id, nsim = 500, nsimra = 500) engen2008I(us, av, nsimra=500) ## S3 method for class 'engenetalI' print(x, ...) ## S3 method for class 'engenetalII' print(x, ...)
us |
a data frame containing the value of numeric habitat variables (columns) in each site (rows) used by the animals. |
av |
a data frame containing the value of numeric habitat variables (columns) in each site (rows) available to the animals. |
id |
a factor with as many elements as there are rows in
|
nsim |
the number of randomizations used in the calculation of the total variance. |
nsimra |
the number of random allocation of ranks used in the calculation of the normal score (see details). |
x |
an object of class |
... |
additional arguments to be passed to other functions (currently unused) |
Engen et al. (2008) proposed an original approach to measure the preference of animals for values of each particular variable of a multivariate set of environmental variables. Their approach was originally developed for the case where there is a sample of used site is for each animal in a sample of identified animals (e.g. using radiotelemetry or GPS), with several sites per animal (i.e., design II according to the classification of Thomas and Taylor, 1990). However, we extended this approach to also include the case where habitat use is described by a sample of used site, with one site per unidentified animal (i.e., design I).
The original approach is the following: first, a normal score
transformation of each habitat variable is performed: for each
variable, the empirical cumulative distribution is computed, by
dividing the rank of the value of each available site by the number of
observations. Note that the ties are ranked randomly. Then, the
inverse of the standard normal integral (see ?qnorm
) of the
cumulative distribution function is computed for all available sites:
this results into a perfectly normal distribution of the habitat
variables for the available sites. Then, the value of the cumulative
distribution – estimated from the available sites – is computed for
each used site. Then, the inverse of the standard normal integral is
computed for each one.
Engen et al. (2008) suppose the following model describing how habitat use results from habitat availability. Let
be the value of a given habitat variable (transformed according to the normal score) for the j-th site used by the i-th animal. Then this value can be described by the model:
where
is the preference for the habitat variable (0 indicates a non-preference),
and
are normal distributions with means equal to zero and variances equal to
and
respectively. Engen et al. give fomula for the estimation of these parameters. Their estimation is done by first estimating the total variance
(this variance is estimated
by sampling randomly one observation per animal – the parameter
nsim
controls the number of samples used in this computation;
see Engen et al. 2008). Note that the correlation between the value
observed for two used units sampled from all the units used by a given
animal is
. A large value of rho indicates a large variation
in the habitat used between animals (or a small within-animal
variation). The main parameter of concern is here the preference.
The function engen2008II
allows to estimate these
parameters.
The function engen2008I
extends this model for design I
studies (a sample of used sites and a sample of available sites,
animals not identified), by considering the following model for these
studies:
where
is the preference for the habitat variable (0 indicates a non-preference), and
are normal distributions with means equal to zero and variances equal to
.
Note that the habitat variables may be correlated on the study area. In this case, observed preference for a given variable may be an artefact of other variables prefered by animals. Supposing that the data.frame containing the
's is a realization of a multivariate normal distribution, we can compute, for each habitat variable and each used site, the *conditional* mean
and *conditional* standard deviation
of this variable at this site, *given* the values of the other habitat variables at this site (this is done using the algorithm described by Ripley, 1987, p.98). We then compute the standardized values
. The preference is then computed using these standardized values.
Because there may be ties in the distribution of values of habitat
variables, the results may vary depending on the random order chosen
for ties when computing the normal scores. Engen et al. recommended
to repeat the function a large number of times, and to use the mean
values as estimates of the parameters. This is what the function
does, and the number of randomization is controlled by the parameter
nsimra
.
Note that all these methods rely on the following hypotheses: (i) independence between animals, (ii) independence between sites, and (iii) all animals are selecting habitat in the same way (in addition to "traditional" hypotheses in these kinds of studies: no territoriality, all animals having equal access to all available resource units, etc., see Manly et al. 2002 for further details).
That the examples below provide an illustration and discussion of interesting and (at first sight) surprising properties of this method.
engen2008I
returns a list of class engenetalI
, and
engen2008II
returns a list of class engenetalII
. Both
types of list contain two elements:
raw |
this is a list containing one data frame per habitat
variable, containing the value of the correlation rho (for
|
results |
a data frame containing the mean values over all the randomizations, of these parameters (columns) for each habitat variable (rows). |
Be patient! these functions can be very long (depending on the number of
sites and on the value of simra
)
Clement Calenge [email protected]
Engen, S., Grotan, V., Halley, D. and Nygard, T. (2008) An efficient multivariate approach for estimating preference when individual observations are dependent. Journal of Animal Ecology, 77: 958–965.
Thomas, D. and Taylor, E. (1990) Study designs and tests for comparing resource use and availability. Journal of Wildlife Management, 54, 322–330.
Ripley, B. (1987) Stochastic Simulation. John Wiley and Sons.
niche
, madifa
,
gnesfa
for another approach to tackle the study of
habitat selection. For categorical variables, see
kselect
## Not run: ################################# ################################# ################################# ## Practical use of engen2008II data(puechabonsp) map <- puechabonsp$map ## Removes the aspect (no factor allowed in the function) slot(map,"data")$Aspect <- NULL ## engen2008II: avail <- slot(map, "data") use <- join(puechabonsp$relocs, map) id <- slot(puechabonsp$relocs, "data")$Name ## This function can be very long: engen2008II(use, avail, id, nsimra=10) ## Practical use of engen2008I data(lynxjura) ma <- lynxjura$map lo <- lynxjura$locs av <- slot(ma, "data") us <- join(lo, ma) us <- us[!is.na(us[,1]),] ## Idem, be patient here: engen2008I(us, av, nsimra=10) ################################# ################################# ################################# ## ## For a deeper discussion on ## this method... a simulation: ################################# ## ## First, simulation of a dataset ## copy and paste this part into R, ## but skip the reading of the ## comments if you are not interested ## into this simulation ## simulate the available points suppressWarnings(RNGversion("3.5.0")) set.seed(235) av <- cbind(rnorm(1000, mean=0, sd=3), rnorm(1000, mean=0, sd=0.5)) tt <- cbind(c(cos(-pi/4),sin(-pi/4)),c(cos(pi/4), sin(pi/4))) av <- as.data.frame(as.matrix(av)%*%tt) ## simulate the used points: we simulate a selection on the first ## principal component of the PCA of the data.frame describing the ## availability. In other words, we simulate the case where the ## habitat selection occurs on the "common part" of the two habitat ## variable (no preference for one particular variable). us <- do.call("rbind", lapply(1:5, function(i) { us1 <- cbind(rnorm(30, mean=rnorm(1, -4, 1), sd=0.5), rnorm(30, mean=rnorm(1, 0, 0.5), sd=0.2)) return(us1%*%tt) })) colnames(us) <- colnames(av) <- c("var1", "var2") id <- gl(5,30) ################################# ## ## Study of the habitat selection on these data ## The data are: ## - us: a matrix containing the used sites for two ## habitat variables ## - av: a matrix containing the available sites for two ## habitat variables ## - id: a vector containing the id of 5 animals ## First illustrate the use and availability of the two variables: plot(av, xlab="Habitat variable 1", ylab="Habitat variable 2", col="grey", pch=16) lius <- split(as.data.frame(us), id) junk <- lapply(1:5, function(i) points(lius[[i]], pch=16, col=i)) ## -----> ***It is very clear that there is a selection***: ## the animals select the low values of both habitat variables. ## (this is what we actually simulated) ## Now perform the method of Engen et al. (2008): engen2008II(us, av, id) ## Surprisingly, the method seems to fail to identify the clear ## habitat selection identified graphically... ## ## In fact, it does not fail: ## this method identifies the part of habitat selection that is clearly ## attributable to a given variable. Here the animals select the ## the common factor expressed in the two variables, and it is impossible ## to identify whether the selection is due only to the variable 1 or to ## the variable 2: it is caused by both variable simultaneously. ## Once the selection on the variable 2 (including the common part) ## has been removed, there is no longer appearant selection on ## variable 1. Once the selection caused by the variable 1 ## (including the common part) has been removed, there is no ## longer selection on variable 2... ## ## For this reason, Engen et al. recommended to use this method ## concurrently with other factor analyses of the habitat selection ## such as madifa, kselect, niche (in ade4 package), etc. ## ## Note also the strong correlation between the value of two random ## points used by a given animal. This indicates a strong variability ## among animals... ## End(Not run)
## Not run: ################################# ################################# ################################# ## Practical use of engen2008II data(puechabonsp) map <- puechabonsp$map ## Removes the aspect (no factor allowed in the function) slot(map,"data")$Aspect <- NULL ## engen2008II: avail <- slot(map, "data") use <- join(puechabonsp$relocs, map) id <- slot(puechabonsp$relocs, "data")$Name ## This function can be very long: engen2008II(use, avail, id, nsimra=10) ## Practical use of engen2008I data(lynxjura) ma <- lynxjura$map lo <- lynxjura$locs av <- slot(ma, "data") us <- join(lo, ma) us <- us[!is.na(us[,1]),] ## Idem, be patient here: engen2008I(us, av, nsimra=10) ################################# ################################# ################################# ## ## For a deeper discussion on ## this method... a simulation: ################################# ## ## First, simulation of a dataset ## copy and paste this part into R, ## but skip the reading of the ## comments if you are not interested ## into this simulation ## simulate the available points suppressWarnings(RNGversion("3.5.0")) set.seed(235) av <- cbind(rnorm(1000, mean=0, sd=3), rnorm(1000, mean=0, sd=0.5)) tt <- cbind(c(cos(-pi/4),sin(-pi/4)),c(cos(pi/4), sin(pi/4))) av <- as.data.frame(as.matrix(av)%*%tt) ## simulate the used points: we simulate a selection on the first ## principal component of the PCA of the data.frame describing the ## availability. In other words, we simulate the case where the ## habitat selection occurs on the "common part" of the two habitat ## variable (no preference for one particular variable). us <- do.call("rbind", lapply(1:5, function(i) { us1 <- cbind(rnorm(30, mean=rnorm(1, -4, 1), sd=0.5), rnorm(30, mean=rnorm(1, 0, 0.5), sd=0.2)) return(us1%*%tt) })) colnames(us) <- colnames(av) <- c("var1", "var2") id <- gl(5,30) ################################# ## ## Study of the habitat selection on these data ## The data are: ## - us: a matrix containing the used sites for two ## habitat variables ## - av: a matrix containing the available sites for two ## habitat variables ## - id: a vector containing the id of 5 animals ## First illustrate the use and availability of the two variables: plot(av, xlab="Habitat variable 1", ylab="Habitat variable 2", col="grey", pch=16) lius <- split(as.data.frame(us), id) junk <- lapply(1:5, function(i) points(lius[[i]], pch=16, col=i)) ## -----> ***It is very clear that there is a selection***: ## the animals select the low values of both habitat variables. ## (this is what we actually simulated) ## Now perform the method of Engen et al. (2008): engen2008II(us, av, id) ## Surprisingly, the method seems to fail to identify the clear ## habitat selection identified graphically... ## ## In fact, it does not fail: ## this method identifies the part of habitat selection that is clearly ## attributable to a given variable. Here the animals select the ## the common factor expressed in the two variables, and it is impossible ## to identify whether the selection is due only to the variable 1 or to ## the variable 2: it is caused by both variable simultaneously. ## Once the selection on the variable 2 (including the common part) ## has been removed, there is no longer appearant selection on ## variable 1. Once the selection caused by the variable 1 ## (including the common part) has been removed, there is no ## longer selection on variable 2... ## ## For this reason, Engen et al. recommended to use this method ## concurrently with other factor analyses of the habitat selection ## such as madifa, kselect, niche (in ade4 package), etc. ## ## Note also the strong correlation between the value of two random ## points used by a given animal. This indicates a strong variability ## among animals... ## End(Not run)
The function gnesfa
allows to perform a general
niche-environment system factor analysis.
gnesfa(dudi, Focus, Reference, centering = c("single", "twice"), scannf = TRUE, nfFirst = 2, nfLast = 0) ## S3 method for class 'gnesfa' print(x, ...)
gnesfa(dudi, Focus, Reference, centering = c("single", "twice"), scannf = TRUE, nfFirst = 2, nfLast = 0) ## S3 method for class 'gnesfa' print(x, ...)
dudi |
an object of class |
Focus |
a vector containing the focus weights |
Reference |
a vector containing the reference weights |
centering |
a character string indicating the type of centering (see details) |
scannf |
a logical value indicating whether the eigenvalues bar plot should be displayed |
nfFirst |
the number of first axes to be kept |
nfLast |
the number of last axes to be kept |
x |
an object of class GNESFA |
... |
further arguments to be passed to other functions |
The GNESFA is an algorithm which generalises several factor analyses of the ecological niche. A table X gives the values of P environmental variables in N resource units (e.g. the pixels of a raster map). A distribution of weights D describes the availability of the resource units to the species (if not specified, these weights are considered to be uniform). Another distribution of weights Dp describes the use of the resource units by the species (for example the proportion of relocations in each pixel of a raster map).
Each environmental variable defines a dimension in a multidimensional space, the ecological space. The N resource units define a cloud of points in this space. Each point is associated to two weights. The GNESFA finds, in the ecological space, the directions on which these two distributions of weights are the most different.
The GNESFA relies on a choice of the analyst, followed by three steps. Before all, the analyst has to choose one distribution of weights as the Reference distribution, and the other one as the Focus distribution; (i) The first table X is centred on the centroid of the Reference distribution; (ii) a principal component analysis of this Reference distribution is performed; (iii) the cloud of points is distorted, so that the Reference distribution takes a standard spherical shape; (iv) a non centred principal component analysis of the Focus distribution allows to identify the directions of the ecological space where the two distributions are the most different.
Depending on the distribution chosen as Reference, this algorithm returns results with different meanings (see examples). This algorithm is closely related to several common analyses of habitat selection/niche (ENFA, MADIFA, Mahalanobis distances, selection ratios, etc.). The examples below give some examples of the mathematical properties of this algorithm.
Note that the function takes a parameter named centering
.
Indeed, two types of centering can be performed prior to the GNESFA.
The choice "single"
consists in the centering of the cloud of
point in the ecological space on the centroid of the Reference
distribution. The choice "twice"
consist to center the cloud
of points on both the centroid of the Reference distribution and the
centroid of the Focus distribution. This is done by projecting the
cloud of points on the hyperplane orthogonal to the marginality vector
(the vector connecting the two centroids. If this choice is done, the
GNESFA is identical to the commonly used Ecological Niche Factor
Analysis (see examples).
gnesfa
returns a list of class gnesfa
containing the
following components:
call |
original call. |
centering |
The type of centering required. |
tab |
a data frame with n rows and p columns. |
Reference |
a vector of length n containing the Reference weights. |
Focus |
a vector of length n containing the Focus weights. |
nfFirst |
the number of kept first axes. |
nfLast |
the number of kept last axes. |
eig |
a vector with all the eigenvalues of the analysis. |
li |
row coordinates, data frame with n rows and nf columns. |
l1 |
row normed coordinates, data frame with n rows and nf columns. |
co |
column scores, data frame with p rows and nf columns. |
cor |
the correlation between the GNESFA axes and the original variables |
Clement Calenge [email protected]
Calenge, C. and Basille, M. (2008) A General Framework for the Statistical Exploration of the Ecological Niche. Journal of Theoretical Biology, 252: 674-685.
Calenge, C., Darmon, G., Basille, M., Loison, A. and Jullien, J.M. (2008) The factorial decomposition of the Mahalanobis distances in habitat selection studies. Ecology, 89, 555-566.
madifa
, mahasuhab
,
enfa
, wi
for closely related methods (see
Examples)
## Not run: ################################################################ ## ## Study of the habitat selection by the chamois in the French ## mountains of Les Bauges ## Loads the data data(bauges) names(bauges) map <- bauges$map locs <- bauges$locs ## displays the data mimage(map) image(map,1) points(locs, pch = 3) ## Prepares the data for the GNESFA: tab <- slot(map, "data") Dp <- slot(count.points(locs,map), "data")[,1] pc <- dudi.pca(tab, scannf = FALSE) ## Example of use with Dp = Reference gn <- gnesfa(pc, Reference = Dp, scannf=FALSE) ## One main axis: barplot(gn$eig) ## The correlation with variables indicate that ## the elevation, the proximity to grass and to ## deciduous forests: s.arrow(gn$cor) ## The factorial map of the niche... scatterniche(gn$li, Dp, pts = TRUE) ## The chamois is rather located at high elevation, ## in the grass, far from deciduous forests ########################################################## ########################################################## ## ## ## Some interesting properties of the GNESFA ## ## ########################################################## ########################################################## ################################ ################################ ## ## Interesting properties of the ## choice: Dp as Reference ## identical to the MADIFA ## (Calenge et al. 2008), ## See the help page of the function madifa ## for other properties) gn <- gnesfa(pc, Reference = Dp, scannf=FALSE, nfFirst = 7) gn ## This is the same as the MADIFA: mad <- madifa(pc, Dp, scannf=FALSE) ## Indeed: plot(gn$li[,1], mad$li[,1]) cor(gn$li[,1], mad$li[,1]) ## And consequently the sum of the squared scores, ## On the axes of the GNESFA... su <- apply(gn$l1,1,function(x) sum(x^2)) ## ... is equal to the Mahalanobis distances between ## the points and the centroid of the niche ## (Clark et al. 1993, see the help page of mahasuhab) su2 <- slot(mahasuhab(map, locs), "data")[,1] ## Indeed: all(su - su2 < 1e-7) plot(su, su2) ################################ ################################ ## ## Centering twice is identical to ## the ENFA (Hirzel et al. 2002, see the help ## page of the function enfa)... ####### ## ## ... If Dp is the Reference: gn <- gnesfa(pc, Reference = Dp, center = "twice", scannf = FALSE) gn enf <- enfa(pc, Dp, scannf = FALSE) plot(enf$li[,2], gn$li[,1]) cor(enf$li[,2], gn$li[,1]) ## The first specialization axis of the ENFA ## is the first axis of the GNESFA! ####### ## ## ... If Dp is the Focus: gn <- gnesfa(pc, Focus = Dp, center = "twice", scannf = FALSE, nfFirst = 6) plot(enf$li[,2], gn$li[,6]) cor(enf$li[,2], gn$li[,6]) ## The first specialization axis of the ENFA ## is the last axis of the GNESFA! ####### ## ## Whatever the distribution chosen as Reference, ## projecting the cloud of points on the hyperplane ## orthogonal to the marginality axis, and performing ## a GNESFA in this subspace is identical to an ENFA! ## The marginality axis of the ENFA is identical ## to the component "projmar" of the GNESFA plot(enf$li[,1],gn$projmar) cor(enf$li[,1],gn$projmar) ################################ ################################ ## ## Interesting properties of the ## case: Dp as Focus, one categorical ## variable. Relationships with the selection ## ratios of Manly et al. (1972, see the ## help page of wi) ## For example, take the Elevation, and ## define a factor with 4 levels elev <- data.frame(el = cut(slot(map, "data")$Elevation, 4)) ## Now, compute the complete disjonctive table dis <- acm.disjonctif(elev) head(dis) ## Now perform the GNESFA with Dp as Focus: pc <- dudi.pca(dis, scannf = FALSE) gn <- gnesfa(pc, Dp, scannf = FALSE, nfFirst = 3) ####### ## ## This analysis is closely related to the concept of ## selection ratios ## Compute the percentage of use of each level: us <- apply(dis, 2, function(x) sum(x*Dp)/sum(Dp)) av <- apply(dis, 2, function(x) sum(x)/length(x)) ## The selection ratios wi <- widesI(us, av)$wi ## Compute the sum of the eigenvalue sum(gn$eig) ## Compute the sum of the selection ratios - 1 sum(wi) - 1 ## In other words, when the GNESFA (Dp as Focus) is ## applied on only one categorical variable, this ## analysis finds a set of axes which partition the ## sum of the selection ratios so that it is maximum ## on the first axes!! ## End(Not run)
## Not run: ################################################################ ## ## Study of the habitat selection by the chamois in the French ## mountains of Les Bauges ## Loads the data data(bauges) names(bauges) map <- bauges$map locs <- bauges$locs ## displays the data mimage(map) image(map,1) points(locs, pch = 3) ## Prepares the data for the GNESFA: tab <- slot(map, "data") Dp <- slot(count.points(locs,map), "data")[,1] pc <- dudi.pca(tab, scannf = FALSE) ## Example of use with Dp = Reference gn <- gnesfa(pc, Reference = Dp, scannf=FALSE) ## One main axis: barplot(gn$eig) ## The correlation with variables indicate that ## the elevation, the proximity to grass and to ## deciduous forests: s.arrow(gn$cor) ## The factorial map of the niche... scatterniche(gn$li, Dp, pts = TRUE) ## The chamois is rather located at high elevation, ## in the grass, far from deciduous forests ########################################################## ########################################################## ## ## ## Some interesting properties of the GNESFA ## ## ########################################################## ########################################################## ################################ ################################ ## ## Interesting properties of the ## choice: Dp as Reference ## identical to the MADIFA ## (Calenge et al. 2008), ## See the help page of the function madifa ## for other properties) gn <- gnesfa(pc, Reference = Dp, scannf=FALSE, nfFirst = 7) gn ## This is the same as the MADIFA: mad <- madifa(pc, Dp, scannf=FALSE) ## Indeed: plot(gn$li[,1], mad$li[,1]) cor(gn$li[,1], mad$li[,1]) ## And consequently the sum of the squared scores, ## On the axes of the GNESFA... su <- apply(gn$l1,1,function(x) sum(x^2)) ## ... is equal to the Mahalanobis distances between ## the points and the centroid of the niche ## (Clark et al. 1993, see the help page of mahasuhab) su2 <- slot(mahasuhab(map, locs), "data")[,1] ## Indeed: all(su - su2 < 1e-7) plot(su, su2) ################################ ################################ ## ## Centering twice is identical to ## the ENFA (Hirzel et al. 2002, see the help ## page of the function enfa)... ####### ## ## ... If Dp is the Reference: gn <- gnesfa(pc, Reference = Dp, center = "twice", scannf = FALSE) gn enf <- enfa(pc, Dp, scannf = FALSE) plot(enf$li[,2], gn$li[,1]) cor(enf$li[,2], gn$li[,1]) ## The first specialization axis of the ENFA ## is the first axis of the GNESFA! ####### ## ## ... If Dp is the Focus: gn <- gnesfa(pc, Focus = Dp, center = "twice", scannf = FALSE, nfFirst = 6) plot(enf$li[,2], gn$li[,6]) cor(enf$li[,2], gn$li[,6]) ## The first specialization axis of the ENFA ## is the last axis of the GNESFA! ####### ## ## Whatever the distribution chosen as Reference, ## projecting the cloud of points on the hyperplane ## orthogonal to the marginality axis, and performing ## a GNESFA in this subspace is identical to an ENFA! ## The marginality axis of the ENFA is identical ## to the component "projmar" of the GNESFA plot(enf$li[,1],gn$projmar) cor(enf$li[,1],gn$projmar) ################################ ################################ ## ## Interesting properties of the ## case: Dp as Focus, one categorical ## variable. Relationships with the selection ## ratios of Manly et al. (1972, see the ## help page of wi) ## For example, take the Elevation, and ## define a factor with 4 levels elev <- data.frame(el = cut(slot(map, "data")$Elevation, 4)) ## Now, compute the complete disjonctive table dis <- acm.disjonctif(elev) head(dis) ## Now perform the GNESFA with Dp as Focus: pc <- dudi.pca(dis, scannf = FALSE) gn <- gnesfa(pc, Dp, scannf = FALSE, nfFirst = 3) ####### ## ## This analysis is closely related to the concept of ## selection ratios ## Compute the percentage of use of each level: us <- apply(dis, 2, function(x) sum(x*Dp)/sum(Dp)) av <- apply(dis, 2, function(x) sum(x)/length(x)) ## The selection ratios wi <- widesI(us, av)$wi ## Compute the sum of the eigenvalue sum(gn$eig) ## Compute the sum of the selection ratios - 1 sum(wi) - 1 ## In other words, when the GNESFA (Dp as Focus) is ## applied on only one categorical variable, this ## analysis finds a set of axes which partition the ## sum of the selection ratios so that it is maximum ## on the first axes!! ## End(Not run)
histniche
draws histograms of the niche-environment system: an
histogram of the available resource units (environment) is drawn on
the same graph as an histogram of the used resource units (i.e. the
niche), for comparison.
histniche(x, pr, type = c("h", "l"), adjust = 1, Acol, Ucol, Aborder, Uborder, Alwd = 1, Ulwd = 1, ylim, ncla = 15, ...)
histniche(x, pr, type = c("h", "l"), adjust = 1, Acol, Ucol, Aborder, Uborder, Alwd = 1, Ulwd = 1, ylim, ncla = 15, ...)
x |
a data frame giving the value of environmental variables (columns) in resource units (rows, e.g. the pixels of a raster map) |
pr |
a vector of integers with the same length as |
type |
what type of plot should be drawn. Possible types are: |
adjust |
if |
Acol |
color for the histograms of the available pixels |
Ucol |
color for the histograms of the used pixels |
Aborder |
if |
Uborder |
if |
Alwd |
if |
Ulwd |
if |
ylim |
the limits for the y axis |
ncla |
The number of classes of the histogram |
... |
further arguments passed to or from other methods |
Mathieu Basille [email protected]
## Not run: data(puechabonsp) cp <- count.points(puechabonsp$relocs, puechabonsp$map) puechabonsp$map histniche(slot(puechabonsp$map, "data"), slot(cp, "data")[,1]) histniche(slot(puechabonsp$map, "data"), slot(cp, "data")[,1], ty="l") ## End(Not run)
## Not run: data(puechabonsp) cp <- count.points(puechabonsp$relocs, puechabonsp$map) puechabonsp$map histniche(slot(puechabonsp$map, "data"), slot(cp, "data")[,1]) histniche(slot(puechabonsp$map, "data"), slot(cp, "data")[,1], ty="l") ## End(Not run)
Performs a multivariate analysis of ecological data (K-select analysis).
kselect(dudi, factor, weight, scannf = TRUE, nf = 2, ewa = FALSE) ## S3 method for class 'kselect' print(x, ...) ## S3 method for class 'kselect' kplot(object, xax = 1, yax = 2, csub = 2, possub = c("topleft", "bottomleft", "bottomright", "topright"), addval = TRUE, cpoint = 1, csize = 1, clegend = 2, ...) ## S3 method for class 'kselect' hist(x, xax = 1, mar=c(0.1,0.1,0.1,0.1), ncell=TRUE, csub=2, possub=c("bottomleft", "topleft", "bottomright", "topright"), ncla=15, ...) ## S3 method for class 'kselect' plot(x, xax = 1, yax = 2, ...) prepksel(sa, hr, locs)
kselect(dudi, factor, weight, scannf = TRUE, nf = 2, ewa = FALSE) ## S3 method for class 'kselect' print(x, ...) ## S3 method for class 'kselect' kplot(object, xax = 1, yax = 2, csub = 2, possub = c("topleft", "bottomleft", "bottomright", "topright"), addval = TRUE, cpoint = 1, csize = 1, clegend = 2, ...) ## S3 method for class 'kselect' hist(x, xax = 1, mar=c(0.1,0.1,0.1,0.1), ncell=TRUE, csub=2, possub=c("bottomleft", "topleft", "bottomright", "topright"), ncla=15, ...) ## S3 method for class 'kselect' plot(x, xax = 1, yax = 2, ...) prepksel(sa, hr, locs)
dudi |
an object of class |
factor |
a factor with the same length as |
weight |
a numeric vector of integer values giving the weight
associated to the rows of |
scannf |
logical. Whether the eigenvalues bar plot should be displayed |
nf |
if |
ewa |
logical. If |
x , object
|
an object of class |
xax |
the column number for the x-axis |
yax |
the column number for the y-axis |
addval |
logical. If |
cpoint |
the size of the points (if 0, the points where no relocations are found are not displayed) |
mar |
the margin parameter (see |
ncell |
logical. If |
csub |
the character size for the legend, used with
|
csize |
the size coefficient for the points |
clegend |
the character size for the legend used by
|
possub |
a character string indicating the sub-title position
|
ncla |
the number of classes of the histograms |
sa |
an object of class |
hr |
an object of class |
locs |
an object of class |
... |
additional arguments to be passed to the generic function
|
The K-select analysis is intended for hindcasting studies of habitat selection by animals using radio-tracking data. Each habitat variable defines one dimension in the ecological space. For each animal, the difference between the vector of average available habitat conditions and the vector of average used conditions defines the marginality vector. Its size is proportional to the importance of habitat selection, and its direction indicates which variables are selected. By performing a non-centered principal component analysis of the table containing the coordinates of the marginality vectors of each animal (row) on the habitat variables (column), the K-select analysis returns a linear combination of habitat variables for which the average marginality is greatest. It is a synthesis of variables which contributes the most to the habitat selection. As with principal component analysis, the biological significance of the factorial axes is deduced from the loading of variables.
prepksel
allows to prepare the data for the kselect analysis
(see examples).
plot.kselect
returns a summary of the analysis: it displays (i)
a graph of the correlations between the principal axes of the PCA of
the objects of class dudi
passed as argument and the factorial
axes of the K-select analysis; (ii) a graph giving the scores of the
habitat variables on the factorial axes of the K-select analysis;
(iii) the barplot of the eigenvalues of the analysis (each eigenvalue
measure the mean marginality explained by the axis; (iv) the
projection of the non-recentred marginality vectors on the factorial
plane (the origin of the arrow indicates the average available habitat
conditions, and the end of the arrow indicates the average used
conditions); (v) the projection of the resource units available to
each animal on the first factorial plane and (vi) the coordinates of
the recentred marginality vectors (i.e. recentred so that they have a
common origin) on the first factorial plane.
kplot.kselect
returns one graph per animal showing the
projections of the available resource units on the factorial plane, as
well as their use by the animal. hist.kselect
does the same
thing, but on one dimension instead of two.
kselect
returns a list of the class kselect
and
dudi
(see dudi
).
Clement Calenge [email protected]
Calenge, C., Dufour, A.B. and Maillard, D. (2005) K-select analysis: a new method to analyse habitat selection in radio-tracking studies. Ecological modelling, 186, 143–153.
s.distri
, and
dudi
for class dudi
.
## Not run: ## Load the data data(puechabonsp) locs <- puechabonsp$relocs map <- puechabonsp$map ## compute the home range of animals (e.g. using the minimum convex ## polygon) pc <- mcp(locs[,"Name"]) ## rasterize it hr <- hr.rast(pc, map) ## Compute the number of relocation in each pixel of the map cp <- count.points(locs[,"Name"], map) ## prepares the data for the kselect analysis x <- prepksel(map, hr, cp) tab <- x$tab ## Example of analysis with two variables: the slope and the elevation. ## Have a look at the use and availability of the two variables ## for the 4 animals tab <- tab[,((names(tab) == "Slope")|(names(tab) == "Elevation"))] tab <- scale(tab) tmp <- split.data.frame(tab, x$factor) wg <- split(x$weight, x$factor) opar <- par(mfrow = n2mfrow(nlevels(x$factor))) for (i in names(tmp)) s.distri(scale(tmp[[i]]), wg[[i]]) par(opar) ## We call a new graphic window x11() ## A K-select analysis acp <- dudi.pca(tab, scannf = FALSE, nf = 2) kn <- kselect(acp, x$factor, x$weight, scannf = FALSE, nf = 2) # use of the generic function scatter scatter(kn) # Displays the first factorial plane kplot(kn) kplot(kn, cellipse = 0, cpoint = 0) kplot(kn, addval = FALSE, cstar = 0) # this factorial plane can be compared with # the other graph to see the rotation proposed by # the analysis graphics.off() # Displays the first factorial axis hist(kn) # Displays the second factorial axis hist(kn, xax = 2) # Summary of the analysis plot(kn) ## End(Not run)
## Not run: ## Load the data data(puechabonsp) locs <- puechabonsp$relocs map <- puechabonsp$map ## compute the home range of animals (e.g. using the minimum convex ## polygon) pc <- mcp(locs[,"Name"]) ## rasterize it hr <- hr.rast(pc, map) ## Compute the number of relocation in each pixel of the map cp <- count.points(locs[,"Name"], map) ## prepares the data for the kselect analysis x <- prepksel(map, hr, cp) tab <- x$tab ## Example of analysis with two variables: the slope and the elevation. ## Have a look at the use and availability of the two variables ## for the 4 animals tab <- tab[,((names(tab) == "Slope")|(names(tab) == "Elevation"))] tab <- scale(tab) tmp <- split.data.frame(tab, x$factor) wg <- split(x$weight, x$factor) opar <- par(mfrow = n2mfrow(nlevels(x$factor))) for (i in names(tmp)) s.distri(scale(tmp[[i]]), wg[[i]]) par(opar) ## We call a new graphic window x11() ## A K-select analysis acp <- dudi.pca(tab, scannf = FALSE, nf = 2) kn <- kselect(acp, x$factor, x$weight, scannf = FALSE, nf = 2) # use of the generic function scatter scatter(kn) # Displays the first factorial plane kplot(kn) kplot(kn, cellipse = 0, cpoint = 0) kplot(kn, addval = FALSE, cstar = 0) # this factorial plane can be compared with # the other graph to see the rotation proposed by # the analysis graphics.off() # Displays the first factorial axis hist(kn) # Displays the second factorial axis hist(kn, xax = 2) # Summary of the analysis plot(kn) ## End(Not run)
The MADIFA allows a factorial decomposition of the Mahalanobis distances. This method is presented here in the framework of niche-environment studies.
predict.madifa
allows the computation of the Mahalanobis
Distances based on a restricted number of factorial axes.
All other functions allow various graphical displays of the results of
the MADIFA.
madifa(dudi, pr, scannf = TRUE, nf = 2) ## S3 method for class 'madifa' print(x, ...) ## S3 method for class 'madifa' scatter(x, xax = 1, yax = 2, pts = FALSE, percent = 95, clabel = 1, side = c("top", "bottom", "none"), Adensity, Udensity, Aangle, Uangle, Aborder, Uborder, Acol, Ucol, Alty, Ulty, Abg, Ubg, Ainch, Uinch, ...) ## S3 method for class 'madifa' hist(x, scores = TRUE, type = c("h", "l"), adjust = 1, Acol, Ucol, Aborder, Uborder, Alwd = 1, Ulwd = 1, ...) ## S3 method for class 'madifa' predict(object, map, nf, ...) s.madifa(x, xax = 1, yax = 2, cgrid = 1, clab = 1, ...) ## S3 method for class 'madifa' plot(x, map, xax = 1, yax = 2, cont = FALSE, ...)
madifa(dudi, pr, scannf = TRUE, nf = 2) ## S3 method for class 'madifa' print(x, ...) ## S3 method for class 'madifa' scatter(x, xax = 1, yax = 2, pts = FALSE, percent = 95, clabel = 1, side = c("top", "bottom", "none"), Adensity, Udensity, Aangle, Uangle, Aborder, Uborder, Acol, Ucol, Alty, Ulty, Abg, Ubg, Ainch, Uinch, ...) ## S3 method for class 'madifa' hist(x, scores = TRUE, type = c("h", "l"), adjust = 1, Acol, Ucol, Aborder, Uborder, Alwd = 1, Ulwd = 1, ...) ## S3 method for class 'madifa' predict(object, map, nf, ...) s.madifa(x, xax = 1, yax = 2, cgrid = 1, clab = 1, ...) ## S3 method for class 'madifa' plot(x, map, xax = 1, yax = 2, cont = FALSE, ...)
dudi |
a duality diagram, an object of class |
pr |
a vector giving the utilization weights associated to each unit |
scannf |
logical. Whether the eigenvalues barplot should be displayed |
nf |
an integer indicating the number of kept factorial axes |
x , object
|
an object of class |
xax |
the column number for the x-axis |
yax |
the column number for the y-axis |
pts |
logical. Whether the points should be drawn. If
|
percent |
100 minus the proportion of outliers to be excluded from the computation of the minimum convex polygons |
clabel |
a character size for the columns |
side |
if |
Adensity |
the density of shading lines, in lines per inch, for
the available pixels polygon. See |
Udensity |
the density of shading lines, in lines per inch, for
the used pixels polygon. See |
Aangle |
the slope of shading lines, given as an angle in degrees (counter-clockwise), for the available pixels polygon |
Uangle |
the slope of shading lines, given as an angle in degrees (counter-clockwise), for the used pixels polygon |
Aborder |
the color for drawing the border of the available pixels
polygon (or of the bars of the histogram) . See
|
Uborder |
the color for drawing the border of the used pixels
polygon (or of the bars of the histogram). See
|
Acol |
the color for filling the available pixels polygon.
if |
Ucol |
the color for filling the used pixels polygon.
if |
Alty |
the line type for the available pixels polygon, as in
|
Ulty |
the line type for the used pixels polygon, as in
|
Abg |
if |
Ubg |
if |
Ainch |
if |
Uinch |
if |
scores |
logical. If |
type |
what type of plot should be drawn. Possible types are: |
adjust |
if |
Alwd |
if |
Ulwd |
if |
cgrid |
a character size, parameter used with par("cex")*
|
clab |
if not NULL, a character size for the labels, used with
|
map |
an object of class |
cont |
logical. Whether contour lines should be added to the maps |
... |
additional arguments to be passed to the functions
|
The Mahalanobis distances are often used in the context of
niche-environment studies (Clark et al. 1993, see the function
mahasuhab
). Each environmental variable defines a dimension in
a multidimensionnal space, namely the ecological space. The
Mahalanobis distance between any resource unit in this space (e.g. the
point defined by the values of environmental variables in a pixel of a
raster map) and the centroid of the niche (the distribution of used
resource units) can be used to give a value of eccentricity to this
point.
For a given distribution of available resource units, for which a
measure of Mahalanobis distances is desired, the MADIFA (MAhalanobis
DIstances Factor Analysis) partitions the ecological space into a set
of axes, so that the first axes maximises the average proportion of
their squared Mahalanobis distances. Note that the sum of the squared
scores of any resource unit on all the axes of the analysis is equal
to the squared Mahalanobis distances for this resource unit. Thus,
the MADIFA partitions the Mahalanobis distances into several axes of
biological meaning (see examples). predict.madifa
allows to
compute approximate Mahalanobis distances from the axes of the
MADIFA.
plot.madifa
returns a graphical summary of the analysis: it
returns graphs of (i) the eigenvalues of the analysis (each eigenvalue
measures the average Mahalanobis distance explained by each factorial
axis); (ii) the scores of the habitat variables (i.e. the coefficients
associated to each environmental variable in the linear combination
defining the axes) - note that as the ecological space is distorted to
"sphericize" the niche, the factorial axes are no longer orthogonals,
and the scores of the variables are distributed within an ellipsoid
instead of an hypersphere of radius equal to one in classical PCA. The
limits of this ellipsoid is displayed on this graph, to see the amount
of distortion done by the analysis (further research needs yet to be
done on this graph); (iii) The projection of the available and used
points on the factorial plane of the MADIFA; (iv) The map of the
Mahalanobis distances computed from the original environmental
variables; (v) the map of the approximated Mahalanobis distances
computed from the two axes displayed in this plot; the correlations
between the original environmental variables and the factorial axes;
(v) the map of the first factorial axis and (vi) the map of the second
factorial axis.
hist.madifa
returns a graph of the niche and the available
resource units on the factorial axes of the analysis.
madifa
returns a list of class madifa
containing the
following components:
call |
original call. |
tab |
a data frame with n rows and p columns. |
pr |
a vector of length n containing the number of points in each pixel of the map. |
nf |
the number of kept factorial axes. |
eig |
a vector with all the eigenvalues of the analysis. |
lw |
row weights, a vector with n components. |
li |
row coordinates, data frame with n rows and nf columns. |
l1 |
row normed coordinates, data frame with n rows and nf columns. |
cw |
column weights, a vector with p components. |
co |
column coordinates, data frame with p rows and nf columns. |
mahasu |
a vector of length n containing the squared Mahalanobis distances for the n units. |
cor |
the correlation between the MADIFA axes and the original variable |
predict.madifa
returns a matrix of class SpatialPixelsDataFrame
.
Clement Calenge [email protected]
Clark, J.D., Dunn, J.E. and Smith, K.G. (1993) A multivariate model of female black bear habitat use for a geographic information system. Journal of Wildlife Management, 57, 519–526.
Calenge, C., Darmon, G., Basille, M., Loison, A. and Jullien J.M. (2008) The factorial decomposition of the Mahalanobis distances in habitat selection studies. Ecology, 89, 555–566.
mahasuhab
for a detailed description of the
Mahalanobis distances, enfa
and gnesfa
for
closely related methods.
## Not run: data(bauges) map <- bauges$map locs <- bauges$loc ## We prepare the data for the MADIFA tab <- slot(map, "data") pr <- slot(count.points(locs, map), "data")[,1] ## We then perform the PCA before the MADIFA pc <- dudi.pca(tab, scannf = FALSE) (mad <- madifa(pc, pr, nf=7, scannf = FALSE)) ######################################### ## ## ## Graphical exploration of the MADIFA ## ## ## ######################################### hist(mad) plot(mad, map) ## this plot represents: ## - the eigenvalues diagram ## - the scores of the columns on the axes ## - a graph of the niche in the available space ## - a map of the Mahalanobis distances computed ## using all environmental variables ## - a map of the Mahalanobis distances computed ## using the two factorial axes used in the ## previous graphs ## - the correlation between habitat variables ## and factorial axes ## - the geographical maps of the two ## factorial axes ## predict with just the first axis pred <- predict(mad, map, nf=1) image(pred) ######################################### ## ## ## Mathematical properties of MADIFA ## ## ## ######################################### ## mad$li is equal to mad$l1, up to a constant (mad$l1 is normed) plot(mad$li[,1],mad$l1[,1]) ## This constant is the square root of the corresponding eigenvalue: ## the variance of mad$l1 is equal to the eigenvalue apply(mad$l1,2,function(x) sum(x^2))/nrow(mad$li) ## the variance of mad$l1 weighted by pr is equal to 1 apply(mad$l1,2,function(x) sum(mad$pr*x^2)/sum(mad$pr)) ## Therefore, the eigenvalues are equal to the average of Mahalanobis ## distance for the available resource units on each axis mean(mahalanobis(matrix(mad$l1[,1], ncol=1), 0, 1)) mad$eig[1] ## Computation of the Mahalanobis distances ma1 <- mahasuhab(map, locs) ## The sum of squared scores for a given Resource unit is equal to the ## Mahalanobis distances ma2 <- apply(mad$l1,1, function(x) sum(x^2)) plot(ma2, slot(ma1, "data")[,1]) ## End(Not run)
## Not run: data(bauges) map <- bauges$map locs <- bauges$loc ## We prepare the data for the MADIFA tab <- slot(map, "data") pr <- slot(count.points(locs, map), "data")[,1] ## We then perform the PCA before the MADIFA pc <- dudi.pca(tab, scannf = FALSE) (mad <- madifa(pc, pr, nf=7, scannf = FALSE)) ######################################### ## ## ## Graphical exploration of the MADIFA ## ## ## ######################################### hist(mad) plot(mad, map) ## this plot represents: ## - the eigenvalues diagram ## - the scores of the columns on the axes ## - a graph of the niche in the available space ## - a map of the Mahalanobis distances computed ## using all environmental variables ## - a map of the Mahalanobis distances computed ## using the two factorial axes used in the ## previous graphs ## - the correlation between habitat variables ## and factorial axes ## - the geographical maps of the two ## factorial axes ## predict with just the first axis pred <- predict(mad, map, nf=1) image(pred) ######################################### ## ## ## Mathematical properties of MADIFA ## ## ## ######################################### ## mad$li is equal to mad$l1, up to a constant (mad$l1 is normed) plot(mad$li[,1],mad$l1[,1]) ## This constant is the square root of the corresponding eigenvalue: ## the variance of mad$l1 is equal to the eigenvalue apply(mad$l1,2,function(x) sum(x^2))/nrow(mad$li) ## the variance of mad$l1 weighted by pr is equal to 1 apply(mad$l1,2,function(x) sum(mad$pr*x^2)/sum(mad$pr)) ## Therefore, the eigenvalues are equal to the average of Mahalanobis ## distance for the available resource units on each axis mean(mahalanobis(matrix(mad$l1[,1], ncol=1), 0, 1)) mad$eig[1] ## Computation of the Mahalanobis distances ma1 <- mahasuhab(map, locs) ## The sum of squared scores for a given Resource unit is equal to the ## Mahalanobis distances ma2 <- apply(mad$l1,1, function(x) sum(x^2)) plot(ma2, slot(ma1, "data")[,1]) ## End(Not run)
This function computes the habitat suitability map of an area for a species, given a set of locations of the species occurences (Clark et al. 1993). This function is to be used in habitat selection studies, when animals are not identified.
mahasuhab(x, pts, type = c("distance", "probability"))
mahasuhab(x, pts, type = c("distance", "probability"))
x |
a raster map of class |
pts |
a data frame with two columns, giving the coordinates of the species locations |
type |
a character string. Whether the raw |
Let assume that a set of locations of the species on an area is available (gathered on transects, or during the monitoring of the population, etc.). If we assume that the probability of detecting an individual is independent from the habitat variables, then we can consider that the habitat found at these sites reflects the habitat use by the animals.
The Mahalanobis distances method has become more and more popular during the past few years to derive habitat suitability maps. The niche of a species is defined as the probability density function of presence of a species in the multidimensionnal space defined by the habitat variables. If this function can be assumed to be multivariate normal, then the mean vector of this distribution corresponds to the optimum for the species.
The function mahasuhab
first computes this mean vector as well
as the variance-covariance matrix of the niche density function, based
on the value of habitat variables in the sample of locations.
Then, the *squared* Mahalanobis distance from this optimum is computed
for each pixel of the map. Thus, the smaller this squared
distance is for a given pixel, and the better is the habitat in this
pixel.
Assuming multivariate normality, squared Mahalanobis distances are
approximately distributed as Chi-square with n degrees of freedom,
where n equals the number of habitat characteristics (see the section
note below on this question). If the
argument type = "probability"
, maps of these p-values are
returned by the function. As such these are the probabilities of a
larger squared Mahalanobis distance than that observed when x is
sampled from the niche.
Returns a raster map of class SpatialPixelsDataFrame
.
The computation of the squared Mahalanobis distances inverts the
variance-covariance matrix of the niche density function (see
?mahalanobis
).
It is therefore important that the habitat variables considered
are not too correlated among each other. When the habitat variables
are too correlated, the variance-covariance matrix is singular and
cannot be inverted.
Note also that it is recommended to scale the variables before the computation, so that they all have the same variance, and therefore the same weight in the analysis (see examples below).
Finally, note that in versions of adehabitatHS prior to 0.3.16,
mahasuhab
incorrectly calculated the probability, when
type = "probability"
, based on a Chi-square with (n-1) degrees
of freedom instead of n degrees of freedom (see Etherington 2019 on
this issue).
Clement Calenge [email protected]
Clark, J.D., Dunn, J.E. and Smith, K.G. (1993) A multivariate model of female black bear habitat use for a geographic information system. Journal of Wildlife Management, 57, 519–526.
Etherington, T.R. (2019) Mahalanobis distances and ecological niche modelling: correcting a chi-squared probability error. Peer J, 7, e6678.
madifa
and dunnfa
for factor
analyses of the Mahalanobis distances, domain
for
another method of habitat suitability mapping,
mahalanobis
for
information on the computation of Mahalanobis distances.
## loads the data data(lynxjura) ma <- lynxjura$map lo <- lynxjura$locs[,1:2] ## We first scale the maps slot(ma, "data") <- dudi.pca(slot(ma, "data"), scannf=FALSE)$tab ## habitat suitability mapping hsm <- mahasuhab(ma, lo, type = "probability") image(hsm) title(main = "Habitat suitability map for the Lynx") points(lo, pch = 3)
## loads the data data(lynxjura) ma <- lynxjura$map lo <- lynxjura$locs[,1:2] ## We first scale the maps slot(ma, "data") <- dudi.pca(slot(ma, "data"), scannf=FALSE)$tab ## habitat suitability mapping hsm <- mahasuhab(ma, lo, type = "probability") image(hsm) title(main = "Habitat suitability map for the Lynx") points(lo, pch = 3)
niche.test
tests for the significance of two parameters of the
ecological niche of a species (marginality and tolerance), using
Monte-Carlo methods. This is a bivariate test.
niche.test(x, pts, nrep = 999, o.include = TRUE, ...)
niche.test(x, pts, nrep = 999, o.include = TRUE, ...)
x |
a raster map of class |
pts |
an object inheriting the class |
nrep |
the number of permutations |
o.include |
logical, passed to |
... |
further arguments passed to |
niche.test
tests the significance of two parameters describing
the ecological niche: the marginality
(squared length of the vector linking the average available habitat
conditions to the average used habitat conditions in the ecological
space defined by the habitat variables), and the tolerance (inertia of
the niche in the ecological space, i.e. the sum over all variables of
the variance of used pixels).
At each step of the randomisation procedure, the test randomly allocates the n points in the pixels of the map. The marginality and the tolerance are then recomputed on this randomised data set.
Actual values are compared to random values with the help of the
function biv.test
.
Returns a list containing the following components:
dfxy |
a data frame with the randomized values of marginality (first column) and tolerance (second column). |
obs |
the actual value of marginality and tolerance. |
biv.test
uses the function kde2d
of the package MASS
.
Mathieu Basille [email protected]
Clement Calenge [email protected]
biv.test
for more details on bivariate tests.
histniche
for the histograms of the variables of the niche.
## Not run: data(chamois) niche <- niche.test(chamois$map, chamois$locs, side = "bottom") names(niche) ## End(Not run)
## Not run: data(chamois) niche <- niche.test(chamois$map, chamois$locs, side = "bottom") names(niche) ## End(Not run)
This data set describes the use and availability of 5 habitat types for 13 pheasants monitored using radio-tracking.
data(pheasant)
data(pheasant)
This list has three components:
studyarea
a data frame giving the proportion of each
habitat type (columns) on the study area. These habitat types are
Scrub
, Broadleaf
, Coniferous
,
Grassland
and Crop
. These proportions are
repeated by rows, for all animals (rows)
mcp
a data frame giving the proportion of each habitat type (columns) in the home range of each animal (rows)
locs
a data frame giving the proportion of relocations
of each animal (rows) reported in 3 of the 5 habitat types
(columns). Indeed, Coniferous
and Crops
were not
used by most of the animals.
Aebischer, N. J., Robertson, P. A. and Kenward, R. E. (1993) Compositional analysis of habitat use from animal radiotracking data. Ecology, 74, 1313–1325.
predict.enfa
computes habitat suitability maps using the
Ecological-Niche Factor Analysis and the Mahalanobis distances
method.
## S3 method for class 'enfa' predict(object, map, nf, ...)
## S3 method for class 'enfa' predict(object, map, nf, ...)
object |
an object of class |
map |
an object of class |
nf |
the number of axes of specialization kept for the
predictions. By default, all axes kept in |
... |
further arguments passed to or from other methods |
The predictions are based on the position of the niche defined by the
ENFA within the multidimensional space of environmental variables. The
ENFA produces row coordinates for each pixel, which are used with the
function mahalanobis
. For each pixel, this function computes the
Mahalanobis distances from the barycentre of the niche.
Actually, the function predict.enfa
is identical to the function
mahasuhab
, except that the habitat suitability map is computed
using the axes of the ENFA, instead of the raw data.
Note that the MADIFA allows a more consistent factorial decomposition of the Mahalanobis distances.
Returns a raster map of class SpatialPixelsDataFrame
.
Mathieu Basille [email protected]
Clark, J.D., Dunn, J.E. and Smith, K.G. (1993) A multivariate model of female black bear habitat use for a geographic information system. Journal of Wildlife Management, 57, 519–526.
Hirzel, A.H., Hausser, J., Chessel, D. & Perrin, N. (2002) Ecological-niche factor analysis: How to compute habitat-suitability maps without absence data? Ecology, 83, 2027–2036.
mahalanobis
for information on the computation of
Mahalanobis distances. mahasuhab
for more details on the
computation of habitat suitability maps using the Mahalanobis distances.
madifa
for a more consistent factorial decomposition of
the Mahalanobis distances
## Not run: data(lynxjura) map <- lynxjura$map ## We keep only "wild" indices. locs <- lynxjura$locs locs <- locs[slot(locs, "data")[,2]!="D",] pr <- slot(count.points(locs, map), "data")[,1] (enfa1 <- enfa(dudi.pca(slot(map, "data"), scannf=FALSE), pr, scannf = FALSE)) ## Compute the prediction pred <- predict(enfa1, map) image(pred) contour(pred, col="green", add=TRUE) points(locs, pch = 3) ## End(Not run)
## Not run: data(lynxjura) map <- lynxjura$map ## We keep only "wild" indices. locs <- lynxjura$locs locs <- locs[slot(locs, "data")[,2]!="D",] pr <- slot(count.points(locs, map), "data")[,1] (enfa1 <- enfa(dudi.pca(slot(map, "data"), scannf=FALSE), pr, scannf = FALSE)) ## Compute the prediction pred <- predict(enfa1, map) image(pred) contour(pred, col="green", add=TRUE) points(locs, pch = 3) ## End(Not run)
This data set stores the results of the monitoring of 6 wild boar at Puechabon (Mediterranean habitat, South of France). These data have been collected by Daniel Maillard (Office national de la chasse et de la faune sauvage).
data(puech)
data(puech)
The list puech
has two components:
puech$relocations
is an object of class
SpatialPixelsDataFrame
containing the relocations of the wild
boar resting sites in summer. It contains the coordinates of the
relocations and the name of the corresponding wild boar.
puech$map
is n object of class SpatialPixelsDataFrame
that describe nine environmental variables on the study area (the
elevation, the tree cover, the shrub cover, the distance to
recreational trails, the distance to crops, the distance to water
points, the grass cover, the slope and the sunshine).
Note that both the maps and the relocations have been slightly "damaged" to preserve the copyright on the data.
Maillard, D. (1996). Occupation et utilisation de la garrigue et du vignoble mediterraneens par le Sanglier. Universite d'Aix-Marseille III: PhD thesis.
This data set contains two data frames describing the use and the availability of 3 elevation classes for 6 wild boars (Sus scrofa L.) monitored using radio-tracking at Puechabon (South of France). These data have been collected by Daniel Maillard (Office national de la chasse et de la faune sauvage).
data(puechdesIII)
data(puechdesIII)
The list puechdesIII
has two components:
The data frame used
describes the number of telemetry
relocations for each of the 6 animals in each of the three elevation
classes.
The data frame available
describes a sample
of random points placed in the areas available to these wild boars
(a buffer area of 200 m around the relocations).
Maillard, D. (1996) Occupation et utilisation de la garrigue et du vignoble mediterraneens par le Sanglier. Universite d'Aix-Marseille III: PhD thesis.
rand.kselect
tests whether the marginality vector of animals
is significantly larger than what is expected under the hypothesis of
random habitat use (third-order habitat selection:
selection by the animals of the relocations within their home range;
the habitat availability is measured for each animal). The effect of
each variable on individual marginality is also tested.
Finally, the pertinence of a K-select analysis is also tested. This is
a randomisation test. The alpha-level of the tests is
ajusted using the Bonferroni inequality.
rand.kselect(dudi, factor, weight, nrep = 200, alpha = 0.05, ewa = FALSE) ## S3 method for class 'rand.kselect' print(x, ...)
rand.kselect(dudi, factor, weight, nrep = 200, alpha = 0.05, ewa = FALSE) ## S3 method for class 'rand.kselect' print(x, ...)
dudi |
an object of class |
factor |
a factor defining the animals identity |
weight |
a weight vector of integer values (number of relocations counted in each resource unit in row of the object dudi) |
nrep |
the number of repetitions of the test |
alpha |
the alpha level for the tests. |
ewa |
logical. If |
x |
an object of the class |
... |
further arguments to be passed to the generic function
|
This test is carried out by simulating a random use of space by
animals. rand.kselect
is closely related to the function
kselect
(same arguments).
At each step of the randomisation procedure, and for each animal, the test randomly allocates the nk relocations (where nk is the sum of the weight vector for the animal k) in the Ik pixels available to this animal (where Ik is the length of the weight vector for animal k).
The length of the marginality vector is recomputed at each step of the randomisation procedure and for each animal. The effect of each variable on the use of pixels by each animal is measured by the criterion "(average habitat variable j used by animal i) minus (average habitat variable j available to animal i)". Finally the value of the first eigenvalue of the K-select analysis provides a criterion to test the pertinence of the K-select analysis.
All these values are then compared to the observed values to assess the significance of theses effects.
Returns an object of class rand.kselect
. This list has three
components:
global |
a vector of length 2 giving the results of the randomisation procedure for the first eigenvalue of the K-select analysis. |
marg |
a matrix giving the significance of the marginality of each animal. |
per.ind |
a list giving the results of the randomisation test for the coordinates of the marginality vector for each animal on each habitat variable. |
alpha |
the alpha level of the tests. |
Clement Calenge [email protected]
Calenge, C., Dufour, A.B. and Maillard, D. (2005) K-select analysis: a new method to analyse habitat selection in radio-tracking studies. Ecological modelling, 186, 143–153.
kselect
to perform a K-select analysis.
## Not run: ## Loads the data data(puechabonsp) locs <- puechabonsp$relocs map <- puechabonsp$map ## compute the home range of animals (e.g. using the minimum convex ## polygon) pc <- mcp(locs[,"Name"]) ## rasterize it hr <- do.call("data.frame", lapply(1:nrow(pc), function(i) { sp::over(map, geometry(pc[i,])) })) names(hr) <- slot(pc, "data")$id coordinates(hr) <- coordinates(map) gridded(hr) <- TRUE ## Compute the number of relocation in each pixel of the map cp <- count.points(locs[,"Name"], map) ## prepares the data for the kselect analysis x <- prepksel(map, hr, cp) tab <- x$tab dud <- dudi.mix(tab, scannf = FALSE, nf = 2) ## the randomisation tests ## be patient, this can be very long on some machines (te <- rand.kselect(dud, x$factor, x$weight, nrep = 500)) ## End(Not run)
## Not run: ## Loads the data data(puechabonsp) locs <- puechabonsp$relocs map <- puechabonsp$map ## compute the home range of animals (e.g. using the minimum convex ## polygon) pc <- mcp(locs[,"Name"]) ## rasterize it hr <- do.call("data.frame", lapply(1:nrow(pc), function(i) { sp::over(map, geometry(pc[i,])) })) names(hr) <- slot(pc, "data")$id coordinates(hr) <- coordinates(map) gridded(hr) <- TRUE ## Compute the number of relocation in each pixel of the map cp <- count.points(locs[,"Name"], map) ## prepares the data for the kselect analysis x <- prepksel(map, hr, cp) tab <- x$tab dud <- dudi.mix(tab, scannf = FALSE, nf = 2) ## the randomisation tests ## be patient, this can be very long on some machines (te <- rand.kselect(dud, x$factor, x$weight, nrep = 500)) ## End(Not run)
randtest.enfa
performs a randomisation test for the Ecological
Niche Factor analysis (ENFA).
## S3 method for class 'enfa' randtest(xtest, nrepet = 999, ...)
## S3 method for class 'enfa' randtest(xtest, nrepet = 999, ...)
xtest |
an object of class |
nrepet |
the number of iterations for the randomisation test |
... |
further arguments to be passed to the generic function
|
This test is carried out by simulating a random distribution of the species occurrences in the pixels of a map.
At each step of the randomisation procedure, the test randomly
allocates the nk occurrences (where nk is the
sum of the occurrence vector pr
of the object of class enfa
)
in the Ik pixels of the focus area (where Ik is the length of this occurrence
vector).
At each step of the procedure, the first eigenvalue of the ENFA performed on the randomised data set is recomputed. This value provides a criterion to test the pertinence of the ENFA analysis.
returns a list of class randtest
Clement Calenge [email protected]
Manly, B.F.J. (1997) Randomization, Bootstrap and Monte Carlo Methods in Biology. London: Chapman & Hall.
Hirzel, A.H., Hausser, J., Chessel, D. and Perrin, N. (2002) Ecological-niche factor analysis: How to compute habitat suitability maps without absence data? Ecology, 83, 2027–2036.
## Not run: data(chamois) locs <- chamois$locs map <- chamois$map ## prepare the data tab <- slot(map, "data") tab$Vegetation <- NULL pr <- slot(count.points(locs, map), "data")[,1] en <- enfa(dudi.pca(tab, scannf=FALSE), pr, scannf = FALSE) (tutu <- randtest(en, nrepet = 100)) plot(tutu) ## End(Not run)
## Not run: data(chamois) locs <- chamois$locs map <- chamois$map ## prepare the data tab <- slot(map, "data") tab$Vegetation <- NULL pr <- slot(count.points(locs, map), "data")[,1] en <- enfa(dudi.pca(tab, scannf=FALSE), pr, scannf = FALSE) (tutu <- randtest(en, nrepet = 100)) plot(tutu) ## End(Not run)
Performs the scatter diagrams of objects of class enfa
.
## S3 method for class 'enfa' scatter(x, xax = 1, yax = 2, pts = FALSE, nc = TRUE, percent = 95, clabel = 1, side = c("top", "bottom", "none"), Adensity, Udensity, Aangle, Uangle, Aborder, Uborder, Acol, Ucol, Alty, Ulty, Abg, Ubg, Ainch, Uinch, ...)
## S3 method for class 'enfa' scatter(x, xax = 1, yax = 2, pts = FALSE, nc = TRUE, percent = 95, clabel = 1, side = c("top", "bottom", "none"), Adensity, Udensity, Aangle, Uangle, Aborder, Uborder, Acol, Ucol, Alty, Ulty, Abg, Ubg, Ainch, Uinch, ...)
x |
an object of class |
xax |
the column number for the x-axis |
yax |
the column number for the y-axis |
pts |
logical. Whether the points should be drawn. If
|
nc |
whether or not the niche center should be displayed |
percent |
100 minus the proportion of outliers to be excluded from the computation of the minimum convex polygons |
clabel |
a character size for the columns |
side |
if |
Adensity |
the density of shading lines, in lines per inch, for the
available pixels polygon. See |
Udensity |
the density of shading lines, in lines per inch, for the
used pixels polygon. See |
Aangle |
the slope of shading lines, given as an angle in degrees (counter-clockwise), for the available pixels polygon |
Uangle |
the slope of shading lines, given as an angle in degrees (counter-clockwise), for the used pixels polygon |
Aborder |
the color for drawing the border of the available pixels
polygon. See |
Uborder |
the color for drawing the border of the used pixels polygon.
See |
Acol |
the color for filling the available pixels polygon.
if |
Ucol |
the color for filling the used pixels polygon.
if |
Alty |
the line type for the available pixels polygon, as in |
Ulty |
the line type for the used pixels polygon, as in |
Abg |
if |
Ubg |
if |
Ainch |
if |
Uinch |
if |
... |
further arguments passed to or from other methods |
scatter.enfa
displays a factorial map of pixels, as well as the
projection of the vectors of the canonical basis multiplied by a
constant of rescaling.
The kept axes for the plot are specified in a corner.
Mathieu Basille [email protected]
Basille, M., Calenge, C., Marboutin, E., Andersen, R. & Gaillard, J.M. (2008) Assessing habitat selection using multivariate statistics: Some refinements of the ecological-niche factor analysis. Ecological Modelling, 211, 233–240.
data(lynxjura) map <- lynxjura$map ## We keep only "wild" indices. locs <- lynxjura$locs locs <- locs[slot(locs, "data")[,2]!="D",] hist(map, type = "l") ## The variable artif is far from symetric ## We perform a square root transformation ## of this variable ## We therefore normalize the variable 'artif' slot(map,"data")[,4] <- sqrt(slot(map,"data")[,4]) hist(map, type = "l") ## We prepare the data for the ENFA tab <- slot(map, "data") pr <- slot(count.points(locs, map), "data")[,1] ## We then perform the PCA before the ENFA pc <- dudi.pca(tab, scannf = FALSE) ## We perform the ENFA (enfa1 <- enfa(pc, pr, scannf = FALSE)) scatter(enfa1)
data(lynxjura) map <- lynxjura$map ## We keep only "wild" indices. locs <- lynxjura$locs locs <- locs[slot(locs, "data")[,2]!="D",] hist(map, type = "l") ## The variable artif is far from symetric ## We perform a square root transformation ## of this variable ## We therefore normalize the variable 'artif' slot(map,"data")[,4] <- sqrt(slot(map,"data")[,4]) hist(map, type = "l") ## We prepare the data for the ENFA tab <- slot(map, "data") pr <- slot(count.points(locs, map), "data")[,1] ## We then perform the PCA before the ENFA pc <- dudi.pca(tab, scannf = FALSE) ## We perform the ENFA (enfa1 <- enfa(pc, pr, scannf = FALSE)) scatter(enfa1)
scatterniche
displays the niche in the Ecological space
(multidimensional space defined by habitat variables).
scatterniche(x, pr, xax = 1, yax = 2, pts = FALSE, percent = 95, clabel = 1, side = c("top", "bottom", "none"), Adensity, Udensity, Aangle, Uangle, Aborder, Uborder, Acol, Ucol, Alty, Ulty, Abg, Ubg, Ainch, Uinch, ...)
scatterniche(x, pr, xax = 1, yax = 2, pts = FALSE, percent = 95, clabel = 1, side = c("top", "bottom", "none"), Adensity, Udensity, Aangle, Uangle, Aborder, Uborder, Acol, Ucol, Alty, Ulty, Abg, Ubg, Ainch, Uinch, ...)
x |
a data frame giving the value of environmental variables (columns) in resource units (rows, e.g. the pixels of a raster map) |
pr |
a vector giving the utilisation weight for each resource unit |
xax |
the column number for the x-axis |
yax |
the column number for the y-axis |
pts |
logical. Whether the points should be drawn. If
|
percent |
100 minus the proportion of outliers to be excluded from the computation of the minimum convex polygons |
clabel |
a character size for the columns |
side |
if |
Adensity |
the density of shading lines, in lines per inch, for the
available pixels polygon. See |
Udensity |
the density of shading lines, in lines per inch, for the
used pixels polygon. See |
Aangle |
the slope of shading lines, given as an angle in degrees (counter-clockwise), for the available pixels polygon |
Uangle |
the slope of shading lines, given as an angle in degrees (counter-clockwise), for the used pixels polygon |
Aborder |
the color to draw the border of the available pixels
polygon. See |
Uborder |
the color to draw the border of the used pixels polygon.
See |
Acol |
the color for filling the available pixels polygon.
if |
Ucol |
the color for filling the used pixels polygon.
if |
Alty |
the line type for the available pixels polygon, as in
|
Ulty |
the line type for the used pixels polygon, as in |
Abg |
if |
Ubg |
if |
Ainch |
if |
Uinch |
if |
... |
further arguments passed to or from other methods |
Mathieu Basille [email protected]
Clement Calenge [email protected]
data(chamois) cpi <- slot(count.points(chamois$locs, chamois$map),"data")[,1] chamois$map tab <- slot(chamois$map, "data") ## we focus on the distance to ecotone and on the slope, ## after centring and scaling (with the help of a PCA) scatterniche(dudi.pca(tab[,2:3], scannf=FALSE)$tab, cpi) scatterniche(dudi.pca(tab[,2:3], scannf=FALSE)$tab, cpi, pts=TRUE)
data(chamois) cpi <- slot(count.points(chamois$locs, chamois$map),"data")[,1] chamois$map tab <- slot(chamois$map, "data") ## we focus on the distance to ecotone and on the slope, ## after centring and scaling (with the help of a PCA) scatterniche(dudi.pca(tab[,2:3], scannf=FALSE)$tab, cpi) scatterniche(dudi.pca(tab[,2:3], scannf=FALSE)$tab, cpi, pts=TRUE)
This data set describes the use and availability of 5 habitat types for
17 squirrels monitored using radio-tracking. See also the dataset
squirreloc
.
data(squirrel)
data(squirrel)
This list has three components:
studyarea
a data frame giving the proportion of each habitat type (columns) on the study area. These proportions are repeated by rows, for all animals
mcp
a data frame giving the proportion of each habitat type (columns) in the home range of each animal (rows)
locs
a data frame giving the proportion of relocations of each animal (rows) reported in each habitat type (columns).
Aebischer, N. J., Robertson, P. A. and Kenward, R. E. (1993) Compositional analysis of habitat use from animal radiotracking data. Ecology, 74, 1313–1325.
This data set contains the trajectories of 15 radio-monitored squirrels, as well as the vector maps of habitat composition.
data(squirreloc)
data(squirreloc)
This data set is a list of two objects:
locs
is a SpatialPointsDataFrame
containing the
relocations of 15 squirrels
map
is an object of class SpatialPolygonsDataFrame
containing the habitat composition of the area. The habitat types
and colour coding are stored in the data.frame
The dataset squirreloc
comes from the Ranges VI software. It
has been used to illustrate the compositional analysis (see
?compana
) and the eigenanalysis of selection ratios (see
?eisera
). See also the dataset squirrel
.
Kenward, R.E., South, A.B. and Walls, S.S. (2003). Ranges6 v1.2 : For the analysis of tracking and location data. Online manual. Anatrack Ltd. Wareham, UK. ISBN 0-9546327-0-2.
This data frame describes the habitat use by the Black Grouse (Tetrao tetrix), the Rock Partridge (Alectoris graeca) and the Rock Ptarmigan (Lagopus mutus), in the Vanoise National Park (French Alps), from 1990 to 2000.
data(vanoise)
data(vanoise)
This data frame has 3110 rows and eight columns describing the habitat
composition
for each occurrence of three species of Galliformes. For each
located occurrence (in rows), the employees of the national park have noted:
the species, the elevation (in metres), the aspect (8 classes), the
habitat type (7 categories, FR
means "fallen rocks") and the
date (season, day, month and year).
Calenge, C., Martinot, J.P. and Lebreton, P. (2003) Ecological niche separation among mountain Galliformes in the Vanoise National Parc. Game and Wildlife Science, 20, 259-285.
These functions compute the resource selection ratios (wi) for design I, II and III data types, with resources defined by several categories. Basic tests are also provided.
widesI(u, a, avknown = TRUE, alpha = 0.05) widesII(u, a, avknown = TRUE, alpha = 0.05) widesIII(u, a, avknown = TRUE, alpha = 0.05) ## S3 method for class 'wiI' print(x, ...) ## S3 method for class 'wiII' print(x, ...) ## S3 method for class 'wiIII' print(x, ...) ## S3 method for class 'wi' plot(x, caxis = 0.7, clab = 1, ylog = FALSE, errbar = c("CI", "SE"), main = "Manly selectivity measure", noorder = TRUE, ...)
widesI(u, a, avknown = TRUE, alpha = 0.05) widesII(u, a, avknown = TRUE, alpha = 0.05) widesIII(u, a, avknown = TRUE, alpha = 0.05) ## S3 method for class 'wiI' print(x, ...) ## S3 method for class 'wiII' print(x, ...) ## S3 method for class 'wiIII' print(x, ...) ## S3 method for class 'wi' plot(x, caxis = 0.7, clab = 1, ylog = FALSE, errbar = c("CI", "SE"), main = "Manly selectivity measure", noorder = TRUE, ...)
u |
for |
a |
for |
avknown |
logical. |
alpha |
the threshold value for the tests and confidence intervals |
x |
an object of class |
caxis |
character size on axes to be passed to
|
clab |
character size of axes labels to be passed to
|
ylog |
logical. If |
errbar |
a character string. Type of error bars: either |
main |
a character string. The title of the graph |
noorder |
logical. If |
... |
additionnal arguments to be passed to the function
|
widesI
may be used to explore resource selection by
animals, when designs I are involved (habitat use and availability are
measured at the population level - animals are not identified). The
function tests habitat selection with the Khi2 of Pearson
and log-likelihood Khi2 (recommended, see Manly et al. 2003). The
Manly selectivity measure (selection ratio = used/available) is
computed, the preference / avoidance is tested for
each habitat, and the differences between selection ratios are computed
and tested.
widesII
computes the selection ratios with design II
data (same availability for all animals, but use measured for each one).
Tests of identical habitat use for all animals, and of habitat
selection are also provided.
widesIII
computes the
selection ratios for design III data (when the use and the
availability are measured for each animal -
see examples on the wild boar below). Habitat selection is tested
using a Chi-square for each animal, and the overall habitat selection
is also tested.
Note that all these methods rely on the following hypotheses: (i) independence between animals, and (ii) all animals are selecting habitat in the same way (in addition to "traditional" hypotheses in these kinds of studies: no territoriality, all animals having equal access to all available resource units, etc., see Manly et al. 2002 for further details).
widesI
returns a list of the class wiI
. widesII
returns a list of class wiII
. widesIII
returns a list of class wiIII
. These objects are all
inheriting from the class wi
. They have the following components:
used.prop |
the proportion of use for each resource type. |
avail.prop |
the proportion of available resources. |
wi |
the Manly selectivity measure (selection ratio: used/available). |
se.wi |
the standard error of the selection ratios. |
comparisons |
a list with the following components:
|
profile |
the profile of preferences: resource types are sorted so that the left type is the most preferred and the right type is the most avoided. Habitats for which the selection ratios are not significant are connected by a line. |
alpha |
the parameter |
avknown |
the parameter |
se.used |
only for designs I, the standard error of the proportion of use. |
se.avail |
only for designs I, the standard error of the available proportion. |
chisquwi |
only for designs I, the results of Chi-Square tests of the hypothesis that the selection ratios are in average equals to zero. |
Bi |
only for designs I, equals to |
Khi2P |
only for designs I, test of random resource use (Pearson statistic). |
Khi2L |
For designs I, test of random resource use (Log-likelihood statistic). For design III, global test of random resource use (Log-likelihood statistic) |
Khi2L1 |
only for designs II, test of identical use of habitat by all animals. |
Khi2L2 |
only for designs II, test of overall habitat selection. |
Khi2L2MinusL1 |
only for designs II, test of hypothesis that animals are on average using resources in proportion to availability, irrespective of whether they are the same or not (= Khi2L2 - Khi2L1). |
wij |
only for designs II and III, a matrix with the selection ratios for all animals and all resource categories. |
ICwiupper |
only for designs II and III, the upper limit of the confidence intervals on the selection ratios. |
ICwilower |
only for designs II and III, the lower limit of the confidence intervals on the selection ratios. |
Khi2Lj |
only for designs III, the test of habitat selection for each animal. |
Clement Calenge [email protected]
Manly B.F.J., McDonald L.L., Thomas, D.L., McDonald, T.L. & Erickson, W.P. (2003) Resource selection by animals - Statistical design and Analysis for field studies. Second edition London: Kluwer academic publishers.
Thomas D. L. and Taylor E. J. (1990) Study designs and tests for comparing resource use and availability. Journal of Wildlife Management, 54, 322–330.
compana
for compositional analysis, and
eisera
to perform an eigenanalysis of selection
ratios.
############################ ## Example of moose (Manly et al., 2003, p.52) ## Known available proportions on design I data moose.avail <- c(0.34, 0.101, 0.104, 0.455) moose.used <- c(25, 22, 30, 40) names(moose.used) <- c("InBurnInterior", "InBurnEdge", "OutOfBurnEdge", "OutOfBurnFurther") names(moose.avail) <- names(moose.used) ## Computation of wi (wiRatio <- widesI(moose.used, moose.avail)) ## plot the values of the selection ratios opar <- par(mfrow=c(2,2)) plot(wiRatio) par(opar) ############################ ## Example of Elk (Manly et al., 2003, p.62) ## Estimated available proportions on design I data elk.avail <- c(15, 61, 84, 40) elk.used <- c(3, 90, 181, 51) names(elk.used) <- c("0%", "1-25%", "26-75%", ">75%") names(elk.avail) <- names(elk.used) ## Computation of wi (wiRatio <- widesI(elk.used, elk.avail, avknown=FALSE)) ## plot the values of the selection ratios opar <- par(mfrow=c(2,2)) plot(wiRatio) par(opar) ############################ ## Example of Bighorn (Manly et al., 2003, p.67) ## Known available proportions on design II data data(bighorn) ## Computation of wi (wi <- widesII(bighorn$used, bighorn$availT, alpha = 0.1)) ## plot the values of the selection ratios opar <- par(mfrow=c(2,2)) plot(wi) ############################ ## Example of Bighorn (Manly et al., 2003, p.74) ## Estimated available proportions on design II data ## Computation of wi (wi <- widesII(bighorn$used, bighorn$availE, avknown = FALSE, alpha = 0.1)) ## plot the values of the selection ratios plot(wi) par(opar) ############################ ## Example of Wild boar ## Estimated available proportions on design III data data(puechdesIII) used <- puechdesIII$used available <- puechdesIII$available ## calculation of the selectio ratios ## with sampled availability (i <- widesIII(used,available, avknown = FALSE, alpha = 0.1)) opar <- par(mfrow = c(2,2)) plot(i) par(opar)
############################ ## Example of moose (Manly et al., 2003, p.52) ## Known available proportions on design I data moose.avail <- c(0.34, 0.101, 0.104, 0.455) moose.used <- c(25, 22, 30, 40) names(moose.used) <- c("InBurnInterior", "InBurnEdge", "OutOfBurnEdge", "OutOfBurnFurther") names(moose.avail) <- names(moose.used) ## Computation of wi (wiRatio <- widesI(moose.used, moose.avail)) ## plot the values of the selection ratios opar <- par(mfrow=c(2,2)) plot(wiRatio) par(opar) ############################ ## Example of Elk (Manly et al., 2003, p.62) ## Estimated available proportions on design I data elk.avail <- c(15, 61, 84, 40) elk.used <- c(3, 90, 181, 51) names(elk.used) <- c("0%", "1-25%", "26-75%", ">75%") names(elk.avail) <- names(elk.used) ## Computation of wi (wiRatio <- widesI(elk.used, elk.avail, avknown=FALSE)) ## plot the values of the selection ratios opar <- par(mfrow=c(2,2)) plot(wiRatio) par(opar) ############################ ## Example of Bighorn (Manly et al., 2003, p.67) ## Known available proportions on design II data data(bighorn) ## Computation of wi (wi <- widesII(bighorn$used, bighorn$availT, alpha = 0.1)) ## plot the values of the selection ratios opar <- par(mfrow=c(2,2)) plot(wi) ############################ ## Example of Bighorn (Manly et al., 2003, p.74) ## Estimated available proportions on design II data ## Computation of wi (wi <- widesII(bighorn$used, bighorn$availE, avknown = FALSE, alpha = 0.1)) ## plot the values of the selection ratios plot(wi) par(opar) ############################ ## Example of Wild boar ## Estimated available proportions on design III data data(puechdesIII) used <- puechdesIII$used available <- puechdesIII$available ## calculation of the selectio ratios ## with sampled availability (i <- widesIII(used,available, avknown = FALSE, alpha = 0.1)) opar <- par(mfrow = c(2,2)) plot(i) par(opar)