TBI  R Documentation 
The function computes and tests Temporal Betadiversity Indices (TBI) between multivariate observations (frequency or presenceabsence data) forming pairs observed at time 1 (T1) and time 2 (T2). The data matrices may contain abundance or presenceabsence data, or other types of frequencylike data (e.g. biomasses). TBI are dissimilarity indices that measure beta differentiation through time. The indices are computed between T1 and T2 for each site. The difference between species (or abundancesperspecies) gains (C/den) and losses (B/den) can be printed out and tested for significance.
TBI(
mat1,
mat2,
method = "%difference",
pa.tr = FALSE,
nperm = 99,
BCD = TRUE,
replace = FALSE,
test.BC = TRUE,
test.t.perm = FALSE,
save.BC = FALSE,
seed. = NULL,
clock = FALSE
)
mat1, mat2 
Two multivariate community composition or gene frequency data matrices (class data.frame or matrix) with the same number of rows and columns. The rows must correspond to the same objects (e.g. sites) and the columns to the same variables, e.g. species or alleles. – The input files are checked for having equal numbers of rows and columns, for rows that are empty in both mat1 and mat2, and for the presence of negative values, which cannot be frequencies. 
method 
One of the following dissimilarity coefficients:

pa.tr 
If 
nperm 
Number of permutations for the tests of significance of the
temporal beta indices and the permutation test of the BC difference. Use

BCD 
If If 
replace 
If 
test.BC 
If 
test.t.perm 
If 
save.BC 
If 
seed. 
Seed for random number generator. If 
clock 
If 
For each object, the function tests the hypothesis (H0) that a species assemblage is not exceptionally different between T1 and T2, compared to assemblages that could have been observed at this site at T1 and T2 under conditions corresponding to H0. If H0 is rejected, the object is recognized as exceptionally different from the other objects for its difference between T1 and T2.
To fix ideas, an example in palaeoecology  A researcher is studying ancient and modern diatom communities in sediment cores. If a site displays an exceptional difference between T1 and T2, the researcher is justified to examine the reason for that difference. It could, for example, be caused by a change in land use at that site, which has caused the difference to be larger than at the other sites, compared to the differences caused by climate change at all sites.
The temporal beta diversity indices available in this function belong to four groups, computed in different ways.
Method
"%difference"
computes the percentage difference index, erroneously
called the BrayCurtis index in some software packages ; it is the
quantitative form of the Sorensen index. Method "ruzicka"
computes the
Ruzicka dissimilarity; this is one of the quantitative coefficients
corresponding to the Jaccard dissimilarity for binary data. When these
indices are used to compute ordinations by principal coordinate analysis, it
is recommended to take the square root of the dissimilarities before the
ordination analysis because these indices do not have the property of being
Euclidean. However, that precaution is not important here; the results of
permutation tests will be the same for these dissimilarities, squarerooted
or not. If pa.tr=TRUE
, either the Sorensen or the Jaccard coefficient
are obtained by computing these two coefficients.
Methods
"chord"
(chord distance), "hellinger"
(Hellinger distance) and
"log.chord"
(log.chord distance) are obtained by transformation of the
species data, as described by Legendre & Borcard (2018), followed by
calculation of the Euclidean distance. For the Hellinger distance, the data
are squarerooted, then subjected to the chord transformation and the
Euclidean distance. For the log.chord distance, the data are transformed by
y' = log(y+1) using function log1p() of R, then subjected to the chord
transformation and the Euclidean distance. These three distances have the
Euclidean property (Legendre & Legendre 2012, Legendre & De Caceres 2013). If
pa.tr=TRUE
, the Ochiai distance for binary data,
sqrt(2)*sqrt(1Ochiai similarity), is obtained from these three coefficients.
Methods "jaccard"
, "sorensen"
, "ochiai"
implement
the Jaccard, Sorensen and Ochiai dissimilarities. For these coefficients, the
data are first transformed to presenceabsence (pa.tr
receives the
value TRUE
), then the dissimilarities are computed using the
corresponding quantitative coefficients (Ruzicka, percentage difference, and
chord).
The Euclidean distance is also available in this function. It is not recommended for community composition or allele frequency data. One can compute it for logtransformed abundance data that do not contain zeros, or very few zeros (short gradients).
The temporal beta indices are tested for significance using permutation tests. The hypotheses are the following:
H0: the site under study (e.g. a species assemblage) is not exceptionally different between T1 and T2, compared to assemblages that could have been observed at this site at T1 and T2 under conditions corresponding to H0. The differences between T1 and T2 all belong to the same statistical population of differences.
H1: the site under study is exceptionally different between times T1 and T2.
In the decomposition of the Ruzicka and percentage difference dissimilarities or their presenceabsence forms (Jaccard, Sorensen), the components B and C are computed as follows:
bj is the part of the abundance of species j that is higher at time 1 than at time 2: bj = (y1j  y2j) if y1j > y2j ; bj = 0 otherwise. B is the sum of the bj values for all species in the group of species under study. It is the unscaled sum of species losses between time 1 and time 2. In the BCD output table BCD.mat, column 1 contains B/den where den is the denominator of the index, i.e. (2A+B+C) for the percentage difference index and (A+B+C) for the Ruzicka index.
cj is the part of the abundance of species j that is higher at time 2 than at time 1: cj = (y2j  y1j) if y2j > y1j ; cj = 0 otherwise. C is the sum of the cj values for all species in the group of species under study. It is the unscaled sum of species gains between time 1 and time 2. In the BCD output table BCD.mat, column 2 contains C/den where den is the denominator of the index, i.e. (2A+B+C) for the percentage difference index and (A+B+C) for the Ruzicka index.
The original values of B and C for each site, without denominator, are also available in the output table BC.
Warning  In real ecological studies, when the TBI test is applied to data where some sites are highly impoverished due to pollution or other extreme environmental situations, this situation may produce sites with very few species (i.e. very low richness) and no species in common for the T1T2 comparisons due to sampling variation at these impoverished sites. The TBI dissimilarity will be high and the test may indicate a significant T1T2 difference if most other sites have higher species richness. This would be a correct statistical outcome for the test. When users of the method identify sites showing significant TBI tests in data, they should check the species richness of these sites at T1 and T2. Interpretation of the test results should be done with caution when high and significant TBI indices are associated with very low richness and no species in common between T1 and T2.
Function TBI returns a list containing the following results:
TBI
The vector of Temporal Betadiversity Indices
(TBI) between observations at times T1 and T2 for each object.
p.TBI
A corresponding vector of pvalues. Significant pvalues
(e.g. p.TBI <= 0.05) indicate exceptional objects for the difference of
their species composition.
p.adj
The pvalues are corrected for multiple testing using
function p.adjust of stats
. The adjustment is done using
method="holm"
, which is the default option of the p.adjust
function.
BCD.mat
An output table with four columns: B/den, C/den,
D=(B+C)/den, and Change. The value den is the denominator of the index,
i.e. (2A+B+C) for the percentage difference index and (A+B+C) for the
Ruzicka index. The decomposition is such that D = B/den + C/den. Columns B
and C indicate which of the D values are associated with large B (losses)
or large C values (gains), before proceeding to the analysis and
interpretation of the D values, using environmental or spatial explanatory
variables, through regression or classification tree analysis. When B > C,
the site has lost species or abundancesperspecies between time 1 and time
2; this is indicated by a "" sign in column Change. On the contrary, if B
< C, the site has gained species or abundancesperspecies between time 1
and time 2; this is indicated by a "+" sign in that column. Sites with
equal amounts of losses and gains are marked with a "0".  The B/den and
C/den values can be plotted in BC plots, which are informative about the
changes that occurred in the data set between the two surveys under study.
 If pa.tr
is TRUE, the B and C components are the numbers of
spepcies losses and gains, and D is either the Sorensen or the Jaccard
dissimilarity.  If BCD=FALSE
, that table is not produced. No table
is (or can be) computed for indices other than the Ruzicka and percentage
difference indices or their binary forms.
BCD.summary
An output table with six columns: mean(B/den);
mean(C/den); mean(D); B/(B+C) (which is mean(B/den) divided by mean(D));
C/(B+C) (which is mean(C/den) divided by mean(D)). These values indicate,
over all sites, which of the processes dominated (loss or gain of species
or abundancesperspecies) when site compositions changed between time 1
and time 2. Change has the same meaning as in table BCD.mat
; the
sign indicates the direction of the mean change over all sites.
t.test_B.C
The results of a paired ttest (parametric) of
significance of the difference between columns C/den and B/den of the
BCD.mat
table. If test.t.perm=TRUE
, the difference between
species gains (C/den) and losses (B/den) is also tested in a permutational
paired ttest and the permutational pvalue is shown in the output table.
This result provides an overall test of the direction of change over all
sites. It helps confirm the asymmetry between species (or
abundancesperspecies) gains (C/den) and species (or
abundancesperspecies) losses (B/den). A star in column p<=0.05 indicates
a significant result of the parametric test at the 0.05 level.
BC
An output table with two columns: B and C. In this table,
the B and C statistics are not divided by a denominator, contrary to the
values B/den and C/den found in the output table BCD.mat
.
Pierre Legendre pierre.legendre@umontreal.ca
Legendre, P. 2019. A temporal betadiversity index to identify exceptional sites in spacetime surveys. Ecology and Evolution (in press).
Legendre, P. & M. De Caceres. 2013. Beta diversity as the variance of community data: dissimilarity coefficients and partitioning. Ecology Letters 16: 951963.
Legendre, P. & D. Borcard. 2018. BoxCoxchord transformations for community composition data prior to beta diversity analysis. Ecography 41: 18201824.
Legendre, P. & L. Legendre. 2012. Numerical Ecology. 3rd English edition. Elsevier Science BV, Amsterdam.
van den Brink, P. J. & C. J. F. ter Braak. 1999. Principal response curves: analysis of timedependent multivariate responses of biological community to stress. Environmental Toxicology and Chemistry 18: 138148.
plot.TBI
if(require("vegan", quietly = TRUE)) {
## Example 1 
## Invertebrate communities subjected to insecticide treatment.
## As an example in their paper on Principal Response Curves (PRC method), van den
## Brink & ter Braak (1999) used observations on the abundances of 178 invertebrate
## species (macroinvertebrates and zooplankton) subjected to treatments in 12 mesocosms by
## the insecticide chlorpyrifos. The mesocosms were sampled at 11 occasions. The data,
## available in the {vegan} package, are logtransformed species abundances, ytranformed =
## log(10*y+1).
## The data of survey #4 will be compared to those of survey #11 in this example.
## Survey #4 was carried out one week after the insecticide treatment, whereas the fauna
## of the mesocosms was considered by the authors to have fully recovered from the
## insecticide treatment at survey #11.
data(pyrifos)
## The mesocosms had originally been attributed at random to the treatments. However,
## to facilitate presentation of the results, they will be listed here in order of
## increased insecticide doses: {0, 0, 0, 0, 0.1, 0.1, 0.9, 0.9, 6, 6, 44, 44} micro g/L.
## Select the 12 data rows of surveys 4 and 11 from the data file and reorder them
ord4 = c(38,39,41,47,37,44,40,46,43,48,42,45)
ord11 = c(122,123,125,131,121,128,124,130,127,132,126,129)
## Run the TBI function
res1 < TBI(pyrifos[ord4,], pyrifos[ord11,], method = "%diff", nperm = 0, test.t.perm = FALSE)
res1$BCD.mat
## Example 2 
## This example uses the mite data available in vegan. Let us pretend that sites 120
## represent T1 and sites 2140 represent T2.
data(mite)
# Run the TBI function
res2 < TBI(mite[1:20,], mite[21:40,], method = "%diff", nperm = 0, test.t.perm = FALSE)
summary(res2)
res2$BCD.mat
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.