practice.Rmd

---
title: "Practice"
output: 
  distill::distill_article:
    toc: true
    toc_depth: 3
    toc_float: true
---

# About


<br />


This tutorial describes the basic workflow showing how to compute step by step
functional diversity (FD) indices in a multidimensional space. It is divided in 
three parts:

1. Computing trait-based distances and the multidimensional functional space 
2. Using the [`mFD`](https://github.com/CmlMagneville/mFD) package to compute 
FD alpha and beta indices and plot them 
([Magneville _et al._ 2021](https://doi.org/10.1111/ecog.05904))
3. Using the [`funrar`](https://github.com/rekyt/funrar) package to compute 
functional rarity indices 
([Violle _et al._ 2017](https://doi.org/10.1016/j.tree.2017.02.002),
[Grenié _et al._ 2017](https://doi.org/10.1111%2Fddi.12629))


<br />


**N.B.** You can chose to do Part 2 and Part 3 in the order that you want but 
Part 1 has to be realized first.


<br />

  
# Prerequisites

Be sure you have followed the [**instructions**](https://frbcesab.github.io/workshop-free/instructions.html)
to set up your system (e.g. `R version >= 3.5`).

If not already done, please install the following R packages: 



```{r echo=TRUE, eval=FALSE}
## CRAN packages ----
pkgs <- c("funrar", "mFD")
install.packages(pkgs)
```


<br /> 
  
  
# Data description
  
The dataset used as study case all along this workshop is the **Fruits dataset** 
based on 25 types of fruits (i.e. species) distributed in 10 fruits baskets 
(i.e. assemblages). Each fruit is characterized by five traits values 
summarized in the following table:
  
<br />
  
| Trait name   | Trait measurement  | Trait type    | Number of classes   | Classes code                       | Unit   |
|:------------:|:------------------:|:-------------:|:-------------------:|:----------------------------------:|:------:|
| Size         | Maximal diameter   | Ordinal       | 5                   | 0-1 ; 1-3 ; 3-5 ; 5-10 ; 10-20     | cm     |
| Plant        | Growth form        | Categorical   | 4                   | tree ; shrub ; vine ; forb         | NA     |
| Climate      | Climatic niche     | Ordinal       | 3                   | temperate ; subtropical ; tropical | NA     |
| Seed         | Seed type          | Ordinal       | 3                   | none ; pip ; pit                   | NA     |
| Sugar        | Sugar              | Continuous    | NA                  | NA                                 | g/kg   |
  
  
<br /> 
  

The use of the `mFD` and `funrar` packages is based on two datasets:

* a **`data.frame`** summarizing traits values for each species called
`fruits_traits` in this tutorial


*Code*
```{r}
## Loading data ----
data("fruits_traits", package = "mFD")

## Removing fuzzy traits in this tutorial ----
fruits_traits <- fruits_traits[ , -c(6:8)]

## Display the table ----
knitr::kable(head(fruits_traits),
             caption = "Species x traits data.frame")
```

<br />

* a **`matrix`** summarizing species assemblages called `baskets_fruits_weights`
in this tutorial. Weights in this matrix can be occurrence data, abundance,
biomass, coverage, etc. The studied example works with biomass (*i.e.* grams of
a fruit in a basket) and this matrix looks as follows:


*Code*
```{r}
## Loading data ----
data("baskets_fruits_weights", package = "mFD")

## Display the table ----
knitr::kable(as.data.frame(baskets_fruits_weights[1:6, 1:6]), 
             centering = TRUE,
             caption = "Species x assemblages matrix based on the **fruits** dataset")
```


# Questions

Using this Practice, we ask the following questions:

* How different are the fruits baskets based on their functional traits?
* __*{{ Question functional rarity }}*__
* __*{{ Question funbiogeo }}*__

To answer these three questions, the first step is to **build a functional space based on species traits** 
on which functional diversity and functional rarity indices will be then computed.

<br/>

# Part 1. Build a functional space using the `mFD` package

<br/>

## 1.0. Compute summaries about your data

<br/>

This part is not developed in this Practice (not enough time to see 
everything ;) ), but it could be useful to know that the `mFD` package can 
compute summaries about your traits or assemblage data. For instance, you can 
compute a matrix of species occurrence in each assemblage (needed in 2.2). 

```{r}
## Summary of the assemblages * species data.frame ----
asb_sp_fruits_summ <- mFD::asb.sp.summary(asb_sp_w = baskets_fruits_weights)
asb_sp_fruits_occ  <- asb_sp_fruits_summ$"asb_sp_occ"

head(asb_sp_fruits_occ)
```

<br/>

## 1.1. What about the traits?

<br/>


The first thing to do before starting analyses is to know your data. To do so,
you must be able to characterize the traits you are using (*i.e.* tell the
package what type of traits you are using). That is why `mFD` package needs a
`data.frame` summarizing the type of each trait (*i.e.* each column of the
`fruits_traits` `data.frame`).

*Code*
```{r}
## Loading data ----
data("fruits_traits_cat", package = "mFD")

## Removing fuzzy traits in this tutorial ----
fruits_traits_cat <- fruits_traits_cat[-c(6:8), ]

## Thus remove the "fuzzy_name" column ----
fruits_traits_cat <- fruits_traits_cat[ , -3]

## Displaying the table ----
knitr::kable(head(fruits_traits_cat), 
             caption = "Traits types based on **fruits & baskets** dataset")
```


The **first column** contains **traits name**. The **second column** contains
**traits type** following this code:

* **N**: nominal trait (factor variable)
* **O**: ordinal traits (ordered variable)
* **C**: circular traits (integer values) (**N.B.** circular traits can not be used in
`mFD` function used to compute functional distance but ok for summary function
and function to group species into Functional Entities)
* **Q**: quantitative traits (numeric values)
* **F**: fuzzy traits (described with several values defined in several columns
in the `fruits_traits` `data.frame`)


<br/>

## 1.2. Computing distances between species based on functional traits

<br/>

The next step toward the computation of functional diversity indices is to
estimate functional traits-based distances between species in order to build the
functional space in which indices will be computed.


To compute trait-based distances, we will use the `mFD::funct.dist()` function
which includes the following arguments:


*Code*
```{r, results = "hide"}
sp_dist_fruits <- mFD::funct.dist(
  sp_tr         = fruits_traits,
  tr_cat        = fruits_traits_cat,
  metric        = "gower",
  scale_euclid  = "scale_center",
  ordinal_var   = "classic",
  weight_type   = "equal",
  stop_if_NA    = TRUE)
```


* `sp_tr` is the species x trait `data.frame`

* `tr_cat` is the `data.frame` summarizing trait type for each trait

* `metric` is a character string referring to the metric used to compute
distances. Two metrics are available and **the choice depends on your traits
data**:

  * if **all traits are continuous** use the **Euclidean distance** (`metric =
  "euclidean"`) and check the [Compute Functional Diversity Indices based on
  Only Continuous
  Traits](https://cmlmagneville.github.io/mFD/articles/Continuous_traits_framework.html)
  tutorial which explains how to build a multidimensional space from traits
  through PCA analysis or considering directly each trait as a dimension.

  * if you have **non-continuous traits** use the **Gower distance** (`metric =
  "gower"`) as this method allows traits weighting. This method can also deal
  with fuzzy traits.

* `scale_euclid` is a character string referring to the way the user wants to
scale **euclidean** traits. You can either chose to scale by range (`range`), use
the center transformation (`center`), use the scale transformation (`scale`), use
the scale-center transformation (`scale_center`) or you can chose not to scale
(`noscale`).

* `ordinal_var` is a character string specifying the method to be used for
ordinal variables (*i.e.* ordered). You can either chose to treat ordinal
variables as continuous variables (with `"classic"` option) or to treat ordinal
variables as ranks (with `metric` or `podani` options, see `mFD::funct.dist()`
help file for detail).

* `weight_type` is a character string referring to the type of method to weight
traits. You can either chose to define weights using the `tr_cat` `data.frame` (cf.
**step 1.1**) (`user` option) or you can chose to give the same weight to all
traits (`equal` option). (**N.B.** Using `mFD`, you can not define weights for fuzzy
traits, use
[`gawdis`](https://cran.r-project.org/package=gawdis) package
instead)

* `stop_if_NA` is a logical value to stop or not the process if the `sp_tr`
`data.frame` contains `NA`. If the `sp_tr` `data.frame` contains `NA` you can either
chose to compute anyway functional distances (but keep in mind that **Functional
measures are sensitive to missing traits!**) or you can delete species with
missing or extrapolate missing traits (see [Johnson _et al._
(2020)](https://onlinelibrary.wiley.com/doi/full/10.1111/geb.13185)).


This function returns a `dist` object with traits-based distances between all
pairs of species:

*Code*
```{r}
## Output of the function mFD::funct.dist() ----
round(sp_dist_fruits, 3)
```


<br/>

## 1.3. Building functional spaces and chosing the best one

<br/>

### 1.3.1. Computing several multimensional functional spaces and assessing their quality

<br/>

In order to generate a multidimensional space in which functional diversity
indices
are computed ([Mouillot _et al._ 2013](https://www.sciencedirect.com/science/article/pii/S0169534712002650), 
we will perform a PCoA using the trait-based distances (and if
required a functional dendrogram).
`mFD` evaluates the quality of PCoA-based multidimensional spaces according
to the deviation between trait-based distances and distances in the functional
space (extension of [Maire _et al._ (2015)](https://onlinelibrary.wiley.com/doi/full/10.1111/geb.12299) framework).
For that, we will use the `mFD::quality.fspaces()` function:

<br />

*Code*
```{r, results = "hide"}
fspaces_quality_fruits <- mFD::quality.fspaces(
  sp_dist             = sp_dist_fruits,
  maxdim_pcoa         = 10,
  deviation_weighting = "absolute",
  fdist_scaling       = FALSE,
  fdendro             = "average")
```


* `sp_dist` is the `dist` object with pairwise trait-based distance between
species as computed in **step 1.2**

* `maxdim_pcoa` is the maximum number of PCoA axes to consider to build
multidimensional spaces. Actually, the maximum number of dimensions considered
depends on the number of PCoA axes with positive eigenvalues.

* `deviation_weighting` refers to the method(s) used to weight the difference
between species pairwise distances in the functional space and trait-based
distances.
**You can chose between**:
  * `absolute`: absolute differences are used to compute the **mean absolute
  deviation (mad)** . It reflects the actual magnitude of errors that will
  affect FD metrics.
  * `squared`: squared differences are used to compute the **root of mean square
  deviation (rmsd)**. This weighting puts more weight to the large deviations
  between trait-based distances and distances in the functional space.
  * Both quality metrics can be used with: `deviation_weighting = c("absolute",
  "squared")`.

* `fdist_scaling` specifies whether distances in the functional space should be
scaled before computing differences with trait-based distances. Scaling ensures
that trait-based distances and distances in the functional space have the same
maximum. Scaling distances implies that the quality of the functional space
accounts for congruence in distances rather than their equality.

**N.B.** The combination of `deviation_weighting` and `fdist_scaling` arguments leads
to **four possible quality metrics**: `mad`, `rmsd`, `mad_scaled` and
`rmsd_scaled`

* `fdendro` specifies the clustering algorithm to compute a functional
dendrogram. `NULL` means no dendrogram computed. The chosen algorithm must be one
of the method recognized by the `stats::hclust()` function from the
[`stats`](https://www.rdocumentation.org/packages/stats) package.


This function returns a list various objects:

<br />

* a `data.frame` gathering for each space (in rows), values of quality metric(s)
(in columns)

*Code*
```{r}
## Quality metrics of functional spaces ----
round(fspaces_quality_fruits$"quality_fspaces", 3)
```

<br />

* `list` with details required for other tasks in **step 1.4** to plot functional
space quality and in **step 1.5** to plot functional space.

<br />

**N.B.** The space with the best quality has the lowest quality metric. Here, thanks
to mad values, we can see that the 4D space is the best one. That is why the
following of this Practice will use this multidimensional space.


<br/>

### 1.3.2. Illustrating the quality of the functional spaces

<br/>

With the `mFD` package, it is possible to illustrate the quality of PCoA-based
multidimensional spaces according to deviation between trait-based distances and
distances in the functional space. For that, we use the
`mFD::quality.fspace.plot()` function with the following arguments:

<br />

*Code*
```{r,  fig.show = 'hide', results = "hide"}
mFD::quality.fspaces.plot(
  fspaces_quality            = fspaces_quality_fruits,
  quality_metric             = "mad",
  fspaces_plot               = c("tree_average", "pcoa_2d", "pcoa_3d", 
                                 "pcoa_4d", "pcoa_5d", "pcoa_6d"),
  name_file                  = NULL,
  range_dist                 = NULL,
  range_dev                  = NULL,
  range_qdev                 = NULL,
  gradient_deviation         = c(neg = "darkblue", nul = "grey80", pos = "darkred"),
  gradient_deviation_quality = c(low = "yellow", high = "red"),
  x_lab                      = "Trait-based distance")
```

* `fspaces_quality` is the output of the `mFD::quality.fspaces()` function
(**step 1.3.1**).

* `quality_metric` refers to the quality metric used. It should be one of the
column name(s) of the table gathering quality metric values (output of
`mFD::quality.fspaces()` called `quality_fspaces`) (here:
`fspaces_quality_fruits$quality_fspaces`) Thus it can be: `mad`, `rmsd`,
`mad_scaled` or `rmsd_scaled` (see **step 1.3.1**)

* `fspaces_plot` refers to the names of spaces for which quality has to be
illustrated (up to 10). Names are those used in the output of
`mFD::quality.fspaces()` function showing the values of the quality metric.

* `name_file` refers to the name of file to save (without extension) if the user
wants to save the figure. If the user only wants the plot to be displayed,
then `name_file = NULL`.

* `range_dist`, `range_dev`, `range_qdev` are arguments to set ranges of panel
axes (check function help for further information).

* `gradient_deviation` and `gradient_deviation_quality` are arguments to set
points colors (check function help for further information).

* `xlab` is a parameter to set x-axis label.

<br />

This function generates a figure with three panels (in rows) for each selected
functional space (in columns). Each column represents a functional space, the
value of the quality metric is written on the top of each column. The x-axis of
all panels represents trait-based distances. The y-axis is different for each
row:

* on the first (top) row, the y-axis represents species functional distances in
the multidimensional space. Thus, the closer species are to the 1:1 line,
the better distances in the functional space fit trait-based ones.
* on the second row, the y-axis shows the raw deviation of species distances in
the functional space compared to trait-based distances. Thus, the raw deviation
reflects the distance to the horizontal line.
* on the third row (bottom), the y-axis shows the absolute or squared deviation of the
("scaled") distance in the functional space. It is the deviation that is taken
into account for computing the quality metric.


*Code*
```{r, fig.height = 7, fig.width = 12, fig.align = "center"}
mFD::quality.fspaces.plot(
  fspaces_quality            = fspaces_quality_fruits,
  quality_metric             = "mad",
  fspaces_plot               = c("tree_average", "pcoa_2d", "pcoa_3d",
                                 "pcoa_4d", "pcoa_5d", "pcoa_6d"),
  name_file                  = NULL,
  range_dist                 = NULL,
  range_dev                  = NULL,
  range_qdev                 = NULL,
  gradient_deviation         = c(neg = "darkblue", nul = "grey80", pos = "darkred"),
  gradient_deviation_quality = c(low = "yellow", high = "red"),
  x_lab                      = "Trait-based distance")
```

<br />


For the 2D space, on the top row there are a lot of points below the 1:1 lines,
meaning that distances are overestimated in this multidimensional space. Looking
at panels, we can see that the 4D space is the one in which points are the
closest to the 1:1 line on the top row,and the closest to the x-axis for the two
bottom rows, which reflects a better quality compared to other functional spaces
/ dendrogram. For the dendrogram, we can see on the top row that species pairs
arrange in horizontal lines, meaning that different trait-based distances have
then the same cophenetic distance on the dendrogram.


<br/>

### 1.3.3. Testing the correlation between functional axes and traits

<br />


`mFD` allows to test for correlations between traits and functional axes and
then illustrate possible correlations (continuous traits = linear model is
computed and r^2^ and associated p-value are returned; non-continuous traits =
Kruskal-Wallis test is computed and eta2 statistic is returned). The function
`mFD::traits.faxes.cor()` allows to test and plot correlation and needs the
following arguments:

* `sp_tr` is the species x traits `data.frame`
* `sp_faxes_coord` is a `matrix` of species coordinates taken from the outputs of
the `mFD::quality.fspaces()` function **with columns representing axes on which
functional space must be computed**. For instance, in this tutorial, we will
plot the functional space for 4 and 10 dimensions (*cf.* the two examples
below). The whole `sp_faxes_coord` can be retrieved through the output of the
`mFD::quality.fspaces()` function:


*Code*
```{r, results = "hide"}
sp_faxes_coord_fruits <- fspaces_quality_fruits$"details_fspaces"$"sp_pc_coord"
```


* `plot` is a logical value indicating whether correlations should be
illustrated or not. If this option is set to `TRUE`, traits-axis relationships
are plotted through scatterplot for continuous traits and boxplot for
non-continuous traits.


The function `mFD::traits.faxes.cor()` works as follows:

<br />

*Code*
```{r, results = "hide"}
fruits_tr_faxes <- mFD::traits.faxes.cor(
  sp_tr          = fruits_traits, 
  sp_faxes_coord = sp_faxes_coord_fruits[ , c("PC1", "PC2", "PC3", "PC4")], 
  plot           = TRUE)
```

We can print only traits with significant effect on position along one of the
axis and look at the plots:


*Code*
```{r, fig.height = 7, fig.width = 12, fig.align = "center"}
## Print traits with significant effect ----
fruits_tr_faxes$"tr_faxes_stat"[which(fruits_tr_faxes$"tr_faxes_stat"$"p.value" < 0.05), ]

## Plot ----
fruits_tr_faxes$"tr_faxes_plot"
```

<br />

We can thus see that **PC1** is mostly driven by *Climate* (temperate on the
left and tropical on the right) and *Plant Type* (forb & shrub on the left vs
tree & vine on the right) and *Size* (large fruits on the right) with weaker
influence of *Seed* (eta2 < 0.25).
Then, **PC2** is mostly driven by *Seed* (no seed on the left and pit seed on
the right) with weaker influence of *Plant Type*. **PC3** is driven by only one
trait, *Size*. And finally **PC4** is mostly driven by *Sugar* (high sugar
content on the right and low sugar content on the left) with a weaker influence
of *Plant Type*.

<br />


## 1.4. Plotting the selected functional space and position of species

<br/>

Once the user has selected the dimensionality of the functional space, `mFD`
allows you to plot the given multidimensional functional space and the position
of species in all 2-dimensions spaces made by pairs of axes.

<br />

The `mFD::funct.space.plot()` function allows to illustrate the position of all
species along pairs of space axes.

<br />

This function allows to plot with many possibilities to change colors/shapes of
each plotted element. Here are listed the main arguments:

* `sp_faxes_coord` is a `matrix` of species coordinates taken from the outputs of
the `mFD::quality.fspaces()` function **with columns representing axes on which
functional space must be computed**. For instance, in this tutorial, we will
plot the functional space for 4 and 10 dimensions (*cf.* the two examples
below). The whole `sp_faxes_coord` can be retrieved through the output of the
`mFD::quality.fspaces()` function:

<br />

*Code*
```{r, results = "hide"}
sp_faxes_coord_fruits <- fspaces_quality_fruits$"details_fspaces"$"sp_pc_coord"
```

<br />


* `faxes` is a `vector` containing names of axes to plot. If set to `NULL`, the
first four functional axes will be plotted.

* `faxes_nm` is a `vector` containing labels of `faxes` (following faxes vector
rank). If `NULL`, labels follow `faxes` vector names.

* `range_faxes` is a `vector` to complete if the user wants to set specific
limits for functional axes. If `range_faxes = c(NA, NA)`, the range is
computed according to the range of values among all axes.

* `plot_ch` is a `logical` value used to draw or not the 2D convex-hull filled
by the global pool of species. Color, fill and opacity of the convex hull
can be chosen through other inputs, please refer to the function's help.

* `plot_sp_nm` is a `vector` containing species names to plot. If `NULL`, no
species names plotted. Name size, color and font can be chosen through other 
inputs, please refer to the function's help.

* `plot_vertices` is a `logical` value used to plot or not vertices with a
different shape than other species. **Be careful** these representations are
2D representations, thus vertices of the convex-hull in the n-multidimensional
space can be close to the center of the hull projected in 2D. Color, fill, 
shape and size of vertices can be chosen through other inputs, please refer to 
the function's help.

* `check_input` is a recurrent argument in the `mFD` package. It defines
whether inputs should be checked before computation or not. Possible error
messages will thus be more understandable for the user than R error messages
(**Recommendation:** set it as `TRUE`).

* other inputs are used to chose color, fill, size, and shape of species from
the global pool, please refer to the function's help.

<br />


Here are the plots for the *fruits & baskets* dataset for the **first four PCoA
axis**:

<br />

*Code*
```{r, fig.height = 7, fig.width = 12, fig.align = "center"}
big_plot <- mFD::funct.space.plot(
  sp_faxes_coord  = sp_faxes_coord_fruits[ , c("PC1", "PC2", "PC3", "PC4")],
  faxes           = c("PC1", "PC2", "PC3", "PC4"),
  name_file       = NULL,
  faxes_nm        = NULL,
  range_faxes     = c(NA, NA),
  plot_ch         = TRUE,
  plot_vertices   = TRUE,
  plot_sp_nm      = NULL,
  check_input     = TRUE)

big_plot$"patchwork"
```

<br />

Here, the convex-hull of the species pool is plotted in white and axis have the
same range to get rid of bias based on different axis scales. Species being
vertices of the 4D convex hull are in purple.

<br/>

# Part 2. Computing and plotting FD indices using the `mFD` package

<br />

The `mFD::alpha.fd.multidim()` function allows computing alpha and beta FD indices. 


<br />

## 2.1. Computing and plotting alpha FD indices


<br />


Using the `alpha.fd.multidim()` function, you can compute up to nine alpha FD indices:

 * `FDis` **Functional Dispersion**: the biomass weighted deviation of species
  traits values from the center of the functional space filled by the assemblage
  *i.e.* the biomass-weighted mean distance to the biomass-weighted mean trait
  values of the assemblage.

  * `FRic` **Functional Richness**: the proportion of functional space filled by
  species of the studied assemblage, *i.e.* the volume inside the convex-hull
  shaping species. To compute `FRic` the number of species must be at least
  higher than the number of functional axis + 1.

  * `FDiv` **Functional Divergence**: the proportion of the biomass supported by
  the species with the most extreme functional traits *i.e.* the ones located
  close to the edge of the convex-hull filled by the assemblage.

  * `FEve` **Functional Evenness**: the regularity of biomass distribution in
  the functional space using the Minimum Spanning Tree linking all species
  present in the assemblage.

  * `FSpe` **Functional Specialization**: the biomass weighted mean distance to
  the mean position of species from the global pool (present in all
  assemblages).

  * `FMPD` **Functional Mean Pairwise Distance**: the mean weighted distance
  between all species pairs.

  * `FNND` **Functional Mean Nearest Neighbour Distance**: the weighted distance
  to the nearest neighbor within the assemblage.

  * `FIde` **Functional Identity**: the mean traits values for the assemblage.
  `FIde` is always computed when `FDis` is computed.

  * `FOri` **Functional Originality**: the weighted mean distance to the nearest
  species from the global species pool.

<br />

**A cheat sheet on alpha FD indices is available [here](biblio/FD_alpha_full.pdf)**

<br />

In this Practice we will only compute five of them.
The function is used as follow:

*Code*
```{r, results = "hide"}
alpha_fd_indices_fruits <- mFD::alpha.fd.multidim(
  sp_faxes_coord   = sp_faxes_coord_fruits[ , c("PC1", "PC2", "PC3", "PC4")],
  asb_sp_w         = baskets_fruits_weights,
  ind_vect         = c("fdis", "fric", "fdiv", 
                       "fspe", "fide"),
  scaling          = TRUE,
  check_input      = TRUE,
  details_returned = TRUE)
```

<br />

The arguments and their use are listed below:

* `sp_faxes_coord` is the species coordinates matrix. This dataframe gathers
only axis of the functional space you have chosen based on **step 4**.

* `asb_sp_w` is the matrix linking species and assemblages they belong to
(summarized in **step 1**).

* `ind_vect` is a vector with names of diversity functional indices to compute.

* `scaling` is a logical value indicating whether indices should be scaled
between 0 and 1. If scaling is to be done, this argument must be set to `TRUE`.

* `check_input` is a recurrent argument in the `mFD` package. It defines whether
inputs should be checked before computation or not. Possible error messages will
thus be more understandable for the user than R error messages
(**Recommendation:** set it as `TRUE`).

* `details_returned` is used if the user wants to store information that are
used in graphical functions. If the user wants to plot FD indices, then
`details_returned` must be set to `TRUE`.

<br />

**N.B.** **Use lowercase letters to enter FD indices names**

<br />

The function has two main outputs:

* a **`data.frame` gathering indices values in each assemblage** (for `FIde`
values, there are as many columns as there are axes to the studied functional
space).

*Code*
```{r}
fd_ind_values_fruits <- alpha_fd_indices_fruits$"functional_diversity_indices"
fd_ind_values_fruits
```

<br />

* a **details list** of `data.frames` and `lists` gathering information such as
coordinates of centroids, distances and identity of the nearest neighbour,
distances to the centroid, etc. The user does not have to directly use it but it
will be useful if FD indices are then plotted. It can be retrieved through:

*Code*
```{r, results = "hide"}
details_list_fruits <- alpha_fd_indices_fruits$"details"
```

<br />

Then, you can plot functional indices using the `mFD::alpha.multidim.plot()` for
up to two assemblages:

<br />

*Code*
```{r, results = "hide",  fig.show = 'hide', message = FALSE, fig.height = 7, fig.width = 12, fig.align = "center"}
plots_alpha <- mFD::alpha.multidim.plot(
  output_alpha_fd_multidim = alpha_fd_indices_fruits,
  plot_asb_nm              = c("basket_1", "basket_5"),
  ind_nm                   = c("fdis", "fric", "fdiv", 
                              "fspe", "fide"),
  faxes                    = NULL,
  faxes_nm                 = NULL,
  range_faxes              = c(NA, NA),
  plot_sp_nm               = NULL,
  save_file                = FALSE,
  check_input              = TRUE) 
```


<br />

This function has a lot of arguments: most of them are graphical
arguments allowing the user to chose colors, shapes, sizes, scales, etc. This
tutorial only presents main arguments. To learn about the use of graphical
arguments, check the function help file. The main arguments of this function are
listed below:

* `output_alpha_fd_multidim` is the output of the ``mFD::alpha.fd.multidim()` 
function.
* `plot_asb_nm` is a vector gathering name(s) of assemblage(s) to plot.
* `ind_vect` is a vector gathering FD indices to plot. Plots are available for
`FDis`, `FIde`, `FEve`, `FRic`, `FDiv`, `FOri`, `FSpe`, and `FNND.`
* `faxes` is a vector containing names of axes to plot. You can only plot from
two to four axes labels for graphical reasons.
* `faxes_nm` is a vector with axes labels if the user ants different axes labels
than `faxes` ones.
* `range_faxes` is a vector with minimum and maximum values for axes.
If `range_faxes = c(NA, NA)`, the range is computed according to the range
of values among all axes, all axes having thus the same range. To have a fair
representation of species positions in all plots, all axes must have the same
range.
* `plot_sp_nm` is a vector containing species names to plot. If `NULL`, then no
name is plotted.
* `check_input` is a recurrent argument in `mFD`. It defines whether inputs
should be checked before computation or not. Possible error messages will thus
be more understandable for the user than R error messages (**Recommendation:**
set it as `TRUE`.
* size, color, fill, and shape arguments for each component of the graphs *i.e.*
species of the global pool, species of the studied assemblage(s), vertices,
centroids and segments. If you have to plot two assemblages, then inputs
should be formatted as follow: `c(pool = ..., asb1 = ..., asb2 = ...)` for
inputs used for global pool and studied assemblages and 
`c(asb1 = ..., asb2 = ...)` for inputs used for studied assemblages only.


<br />

Then, using these arguments, here are the output plots for the **fruits &
baskets** dataset:

<br />

* `FRic` representation: the colored shapes reflect the convex-hull of the 
studied assemblages and the white shape reflects the convex-hull of the global 
pool of species:

*Code*
```{r, fig.height = 15, fig.width = 20, fig.align = "center", warning = FALSE}
plots_alpha$"fric"$"patchwork"
```


* `FDiv` representation: the gravity centers of **vertices** (i.e. species with the
most extreme functional traits) of each assemblages are plotted as a square and 
a triangle. The two colored circles represent the mean distance of species to 
the gravity center for each assemblage. Species of each assemblage have 
different size given their relative weight into the assemblage.

*Code*
```{r, fig.height = 15, fig.width = 20, fig.align = "center", warning = FALSE}
plots_alpha$"fdiv"$"patchwork"
```


* `FSpe` representation: colored traits represent distances of each species from
a given assemblage to the center of gravity of the global pool (i.e center of 
the functional space). the center of gravity is plotted with a purple diamond.
Species of each assemblage have different size given their relative weight into 
the assemblage.

*Code*
```{r, fig.height = 15, fig.width = 20, fig.align = "center", warning = FALSE}
plots_alpha$"fspe"$"patchwork"
```


* `FDis` representation: colored traits represent distances of each species from
a given assemblage to the center of gravity of species of the assemblage
(defined by FIde values). The center of gravity of each assemblage is plotted 
using a square and a triangle. Species of each assemblage havedifferent size 
given their relative weight into the assemblage.

*Code*
```{r, fig.height = 15, fig.width = 20, fig.align = "center", warning = FALSE}
plots_alpha$"fdis"$"patchwork"
```


* `FIde` representation:colored lines refer to the weighted average position of 
species of each assemblage along each axis. Species of each assemblage have 
different size given their relative weight into the assemblage.

*Code*
```{r, fig.height = 15, fig.width = 20, fig.align = "center", warning = FALSE}
plots_alpha$"fide"$"patchwork"
```


**N.B.** Using the `mFD` package, you can plot more than two assemblages but not 
with the `alpha.multidim.plot()` function. There are several specific functions 
for each step of the plot: build the background of the plot (`background.plot()`),
plot the pool of species you are working on (`pool.plot()`), 
plot species from the studied assemblages (`species.plot()`) function and lastly 
plot the wanted metric using related function (`fric.plot()`, `fdiv.plot()`, 
`fide.plot()`,`fdis.plot()`, `feve.plot()`, `fnnd.plot()`, `fori.plot()`, 
`fspe.plot()`). Plots for different axes combination can be gathered into a 
single plot using the `panels.to.patchwork()` function.


<br/>

## 2.2. Computing and plotting beta FD indices

<br/>

**N.B.** Some Mac OS X 10.15 may encounter some issues with the beta_*() 
functions.


`mFD` package allows you to compute beta diversity indices for each assemblage
pairs following [Villeger _et al._
2013](https://onlinelibrary.wiley.com/doi/full/10.1111/geb.12021). For that we
will use the `mFD::beta.fd.multidim()` function. This function can compute two 
families of functional beta diversity indices, either **Jaccard** or **Sorensen**.


In this example, we will use Jaccard index. For each assemblages pair, the
dissimilarity index is decomposed into two additive components: **turnover** and
**nestedness-resultant**. 

**NB** The **turnover** component is the highest if there
is no shared traits combination between the two assemblages.
The **nestedness** component is the highest if one assemblage hosts a small
subset of the functional strategies present in the other.


The `mFD::beta.fd.multidim()` function has the main following arguments:

*Code*
```{r, results = "hide", message=FALSE}
beta_fd_indices_fruits <- mFD::beta.fd.multidim(
      sp_faxes_coord   = sp_faxes_coord_fruits[ , c("PC1", "PC2", "PC3", "PC4")],
      asb_sp_occ       = asb_sp_fruits_occ,
      check_input      = TRUE,
      beta_family      = c("Jaccard"),
      details_returned = TRUE)
```


* `sp_faxes_coord` is the species coordinates matrix. This dataframe gathers
**only** axis of the functional space you have chosen based on **step 4**.

* `asb_sp_occ` is the matrix of occurrence (coded as 0/1) of species assemblages
(summarized in **step 1**).

* `check_input` is a recurrent argument in the `mFD` package. It defines whether
inputs should be checked before computation or not. Possible error messages will
thus be more understandable for the user than R error messages (**Recommendation:**
set it as `TRUE`.

* `beta_family` a character string for the type of beta-diversity index to
compute, it can either be `Jaccard` or `Sorensen`.

* `details_returned` is a logical value indicating whether details of outputs
must be stored. It should be stored if you plan to use the graphical function to
illustrate beta diversity indices thereafter.

* There are also other arguments for parallelisation options. Check the function
help file for more explanation.

<br />

The function returns a list containing:

* a dist object with beta indices values for each pair of assemblages:

*Code*
```{r}
head(beta_fd_indices_fruits$"pairasb_fbd_indices", 10)
```


* a list containing details such as inputs, vertices of the global pool and of
each assemblage and FRic values for each assemblage

*Code*
```{r}
beta_fd_indices_fruits$"details"
```


* a vector containing the `FRic` value for each assemblage retrieved through
the `details_beta` list:

*Code*
```{r}
beta_fd_indices_fruits$"details"$"asb_FRic"
```

<br />

* a list of vectors containing names of species being vertices of the convex
hull for each assemblage retrieved through the `details_beta` list:

```{r}
beta_fd_indices_fruits$"details"$"asb_vertices"
```


Then, the package allows the user to illustrate functional beta-diversity
indices for a pair of assemblages in a multidimensional space using the
`mFD::beta.multidim.plot()` function. The output of this function is a figure
showing the overlap between convex hulls shaping each of the two species
assemblages.


The plotting function has a large number of arguments, allowing the user to
chose graphical options. Only main arguments are listed below:

*Code*
```{r, results = "hide", warning = FALSE, fig.height = 7, fig.width = 12, fig.align = "center"}
beta_plot_fruits <- mFD::beta.multidim.plot(
  output_beta_fd_multidim = beta_fd_indices_fruits,
  plot_asb_nm             = c("basket_1", "basket_4"),
  beta_family             = c("Jaccard"),
  plot_sp_nm              = c("apple", "lemon", "pear"),
  faxes                   = paste0("PC", 1:4),
  name_file               = NULL,
  faxes_nm                = NULL,
  range_faxes             = c(NA, NA),
  check_input             = TRUE)
```


* `output_beta_fd_multidim` is the output of the `mFD::beta.fd.multidim()`
function retrieved before as `beta_fd_indices`.

* `plot_asb_nm` is a vector containing the name of the two assemblages to plot.
Here plots of indices will be shown for *basket_1* and *basket_4*.

* `beta_family` refers to the family of the plotted index. It must be the same
as the family chosen to compute beta functional indices values with the
`mFD::beta.fd.multidim()` function.

* `plot_sp_nm` is a vector containing the names of species the user want to
plot, if any. If no the user does not want to plot any species name, then this
argument must be set up to `NULL`. Here, *apple*, *cherry* and *lemon* will be
plotted on the graph.

* `faxes` is a vector containing the names of the functional axes of the plotted
functional space. Here, the figure will be plotted for *PC1*, *PC2* and *PC3*.
This function allows you to plot between two and four axes for graphical
reasons.

* `name_file` is a character string with the name of the file to save the figure
(without extension). If the user does not want to save the file and only display
it, this argument must be set up to `NULL`.

* `faxes_nm` is a vector containing the axes labels for the figure if the user
wants to set up different labels than those contained in `faxes`.

* `range_faxes` is a vector with minimum and maximum values of functional axes.
To have a fair representation of the position of species in all plots, axes
should have the same range. If the user wants the range to be computed according
to the range of values among all axes, this argument must be set up to
`c(NA, NA)`.

* `check_input` is a recurrent argument in the `mFD` package. It defines whether
inputs should be checked before computation or not. Possible error messages will
thus be more understandable for the user than R error messages
(**Recommendation:** set it as `TRUE`)

* Others arguments to set up colors, shapes, sizes and, text fonts are also
available. For more information about them, read the function help file.


Then, the function returns each graph for each functional axes combination and
also a multipanel plot with all combinations of axes and the graph caption.
Here is the multipanel for the fruits exaample:


```{r, fig.height = 16, fig.width = 16, fig.align = "center", warning = FALSE}
beta_plot_fruits$"patchwork"
```

For each assemblage, the associated convex hull is plotted in a different colour
and indices values are printed on the right corner of the plot. Vertices of the
convex hull of a given assemblage can be plotted with a different symbol such as
in this example. Species of all assemblages are plotted with gray cross and the
associated convex hull is plotted in white.

<br/>

# Part 3. Functional rarity

In this part we'll be using our built functional space to compute functional
originality indices. The questions we'll answer are:

* What are functionally distinct fruits compared to all fruits together?
* Are some of these fruits distinct at regional scale but not when looking at specific baskets?

## 3.1 Components of functional rarity: trait originality and rarity

Functional originality indices (also named functional rarity indices) measures
how different are the trait of a target species compared to other ones among
a set. Their idea is analogous to the rarity concept 
in terms of abundance: some species are rare because they show low abundance,
some are common because of their high abundance. Well, functionally original
species are original because they display trait values distant from most other
species.

The indices were initially proposed as complimentary facets: one for "classical"
rarity and one for trait "rarity" (= originality) [Violle _et al._ 2017](https://doi.org/10.1016/j.tree.2017.02.002). With this framework,
a species could be either rare or common, with either original or common traits.
All combinations are possible. For example a species can be rare in abundance
but have common traits. A species could also be common in abundance but original
in terms of traits.

Violle _et al._ ([2014](https://doi.org/10.1016/j.tree.2017.02.002)) further
declined the framework along two spatial scales, local and regional. The local
scale is scale of the site or the assemblage, where species interact directly.
The regional scale is the species pools scale, containing all species occurring
in a given set of assemblages.

This initial framework thus contained 16 different possibilities between
the different kinds of rarity (trait or abundance) and the spatial scales.


## 3.2 Different indices of functional rarity

One specificity of functional rarity indices is that we compute them
on a species basis and not per sites. This means that for each species we will
get one value.

### 3.2.1 Functional originality indices

The two proposed indices of functional originality are
**functional distinctiveness** (*Di*) and **functional uniqueness** (*Ui*).
They were initially envisioned at different spatial scales, with functional
distinctiveness computed at local scale (and considering abundance) while
functional uniqueness was to be computed at regional scales (without considering
abundances).

![`r fontawesome::fa_i("caret-right")` &nbsp;Scheme showing how Distinctiveness (Di) and Uniqueness (Ui) are computed](img/di_ui_scheme.jpg)

**Functional distinctiveness** is the mean of dissimilarity of 
the focal species to all the other species of the set of interest. It can be
abundance-weighted if needed.

**Functional uniqueness** is the smallest dissimilarity that exists between
the focal species and the all other species in the set. It does not consider
the abundance of any species.

### 3.2.1 Rarity indices indices

The two proposed indices of rarity are **scarcity** (*Si*) and
**restrictedness**. They are envisioned at different spatial scales because
scarcity is directly computed from local species relative abundance while
restrictedness is computed at regional scales considering all sites a species
can occupy.

**Scarcity** is proportional to the relative abundance of the species. It gets
close to one when the species is (relatively) rare and close to 0 when its
dominant.

**Restrictedness** is 1 minus the ratio of sites a species occupy over the total
number of sites.

## 3.3 Computing functional rarity

To compute functional rarity indices we will use two object: a functional dissimilarity matrix and a site-species matrix. We will be reusing the ones
we defined in part 1 for that.

Because the `funrar` package follows strictly the framework defined by Violle _et al._ (2017), we have to consider the spatial scale at which we compute distinctiveness and uniqueness.

### 3.3.1 Functional originality at regional scale

At regional scale, we can compute distinctiveness using the `distinctiveness_global()` function. The first argument should be a functional
dissimilarity matrix and the second argument gives
the name of the output column.

```{r di-regional}
library("funrar")

sp_di <- distinctiveness_global(sp_dist_fruits, di_name = "distinctiveness")

head(sp_di)
dim(sp_di)
```

We get one value of distinctiveness per species. It only considers
the functional dissimilarity of all species in the dissimilarity matrix without
considering their spatial distributions. We get one value of distinctiveness per
species.

```{r di-regional-summary, results='hide'}
summary(sp_di)

quantile(sp_di$distinctiveness, probs = seq(0, 1, by = 0.1))

subset(sp_di, distinctiveness >= 0.48)
```
 
When looking at the most distinct fruits we see that banana, litchi,
and pineapple, are the most distinct from their characteristics.

For the choice or dissimilarity matrix we can use the raw dissimilarity
matrix computed directly on raw traits values among species, as we did here.
Another option would be to compute a new functional dissimilarity matrix based
on the selected functional axes. One advantage of the latter is that it already
takes into account the correlation between traits.

Let's recompute regional functional distinctiveness based on the four selected
functional axes. Because the space comes from a PCA, we can directly use
euclidean distance.

```{r di-global-alt}
new_dissim <- dist(sp_faxes_coord_fruits[, c("PC1", "PC2", "PC3", "PC4")])

sp_di_alt <- distinctiveness_global(new_dissim, di_name = "alt_di")
```

We can now compare both distinctiveness values.

```{r di-comparison}
sp_all_di <- merge(sp_di, sp_di_alt, by = "species")

plot(sp_all_di$distinctiveness, sp_all_di$alt_di)
cor.test(sp_all_di$distinctiveness, sp_all_di$alt_di)
```

Both seems very correlated, so in our case using either one should be fine.
However, it can be better to use dissimilarity based on a reduced number of
well-defined axes because: (1) there are more interpretable thanks to the multivariate analysis, (2) the first one contain de most information,
(3) they explicitly take into account potentially strong correlations between
provided traits.
We'll stick here with raw dissimilarity for the sake of simplicity.

To compute uniqueness at regional scale we also need the regional level
functional dissimilarity matrix with the `uniqueness()` function, and
the site-species matrix:

```{r ui-regional, results='hide'}
sp_ui <- uniqueness(
  pres_matrix = baskets_fruits_weights,
  as.matrix(sp_dist_fruits)
)

head(sp_ui)
quantile(sp_ui$Ui, probs = seq(0, 1, by = 0.1))

subset(sp_ui, Ui >= 0.21)
```

Note that we have to transform the dissimilarity object `sp_dist_fruits` into
a matrix explicitly so that the function `uniqueness()` works.

Based on these results we see that banana, passion fruit, and pineapple are
the most isolated fruits in the functional space. Meaning that they have
the most distant nearest neighbors.

So we see that some species can be regionally distinct and regionally unique
(banana and pineapple), while some can be not so regionally distinct but
regionally unique (passion fruit).


### 3.3.2 Functional originality at local scale

At local scale we can use similar logic but we will obtain one value per species
per site.

For distinctiveness we can use the `distinctiveness()` function:

```{r di-local}
sp_local_di <- distinctiveness(
  baskets_fruits_weights, as.matrix(sp_dist_fruits)
)

sp_local_di[1:6, 1:6]

identical(dim(sp_local_di), dim(baskets_fruits_weights))
```
**Note**: because we provided a site-species matrix that contained "abundance"
values (here the weights of the fruits in our baskets), the `distinctiveness()`
is going to use them in its computation. If we want it not to consider
abundances we need to convert the site-species matrix into a presence-absence
matrix.

We get exactly the same number of sites and species as provided in the
lists of baskets. Looking at the first column we see that the distinctiveness
of apple varies a lot. Let's see if can see how much it varies and why.

```{r di-local-apple}
sp_local_di[, 1]
```
It varies from around 0.11 in basket 4 to up to 0.41 in basket 6. We can lookup
the composition of these baskets to understand why:

```{r di-local-basket-compo}
baskets_fruits_weights[c(4, 6),]
```
We can see that basket 4 contains fruits of similar size, origin and sweetness of apples: peach, pear, plums. This decreases the local distinctiveness of apples.
While basket 6 contains more tropical fruits that are quite different in terms
of size and sugar from apples.

To compute uniqueness at the site scale, we must use a more complex expression
as it was not envisioned for local computation:

```{r ui-local}
basket_ui <- apply(
  baskets_fruits_weights, 1,
  function(single_site, dist_m) {
    single_site = single_site[single_site > 0 & !is.na(single_site)]
    uniqueness(t(as.matrix(single_site)), dist_m)
  }, dist_m = as.matrix(sp_dist_fruits)
)

head(basket_ui[1])
```
As we had to manually build the function to compute the local uniqueness
the results are strangely formatted.

We provide here a function that can help them to be more easily read:

```{r ui-local-transform}
basket_ui <- lapply(names(basket_ui), function(x) {
  single_basket = basket_ui[[x]]
  single_basket$site = x
  
  return(single_basket)
})

basket_ui <- do.call(rbind, basket_ui)[, c(3, 1, 2)]
```


Then we can again look at the apple to see how its uniqueness varies across
baskets.

```{r ui-local-apple}
subset(basket_ui, species == "apple")
```

We see that apple have the same uniqueness everywhere but in baskets 6 and 7.
This means that it has the same nearest neighbor (in functional space) across
all of these baskets. But that its closest neighbor is further away in baskets
6 and 7.


### 3.3.3 Rarity indices

Because rarity indices included in this framework are directly tied to specific
spatial scales, we will compute them only at their respective spatial scales.

Similarly to functional originality indices, the functions to compute rarity
indices follow their names. So we can use the `scarcity()` function to compute
scarcity:

```{r scarcity}
si = scarcity(baskets_fruits_weights)
```
Because scarcity needs relative abundances we get prompted by a message asking
us to check if our data is relative. We will use `make_relative()` to transform
the site-species matrix to a relative abundance matrix.

```{r scarcity-working}
rel_weights = make_relative(baskets_fruits_weights)

si = scarcity(rel_weights)

si[1:4, 1:4]
summary(si)
```

We obtain an object in a similar shape as the site-species matrix with local
scarcity values for each species at each site. It gives `NA` for species that
are absent.

To compute restrictedness we can use the `restrictedness()` function.

```{r restrictedness}
ri = restrictedness(baskets_fruits_weights)

head(ri)

summary(ri)
```

This time, because it's a regional index it outputs one value per species as
shown in the `Ri` column.

## 3.4 Plotting functional rarity

`funrar` does not come with visualization functions so that you can use whatever
tool you prefer to plot the metrics. Here we are going to use `ggplot2`
visualize the results. However, feel free to use anything plotting packages
you would like. If you're not fluent in `ggplot2` don't worry this section will
mostly be about interpreting the figures rather than commenting on code.

### 3.4.1 Plotting functional originality

We already saw in part 1 how to visualize the functional space using
the pre-made functions of `mFD`. Here we will use our own functions to be able
to color species in function of their functional originality.

```{r funrar-summary}
# Make a summary data.frame
sp_coord_di_ui <- as.data.frame(sp_faxes_coord_fruits[, 1:2])
sp_coord_di_ui$species <- rownames(sp_coord_di_ui)
rownames(sp_coord_di_ui) <- NULL
sp_coord_di_ui <- sp_coord_di_ui[, c(3, 1, 2)]
sp_coord_di_ui <- merge(sp_coord_di_ui, sp_di, by = "species")
sp_coord_di_ui <- merge(sp_coord_di_ui, sp_ui, by = "species")
```


```{r funrar-plot-di-ui-regional}
library("ggplot2")

plot_reg_distinctiveness <- ggplot(sp_coord_di_ui, aes(PC1, PC2)) +
  geom_hline(yintercept = 0, linetype = 2) +
  geom_vline(xintercept = 0, linetype = 2) +
  geom_point(aes(color = distinctiveness)) +
  ggrepel::geom_text_repel(aes(label = species)) +
  scale_color_viridis_c("Functional\nDistinctiveness") +
  theme_bw() +
  theme(aspect.ratio = 1)

plot_reg_uniqueness <- ggplot(sp_coord_di_ui, aes(PC1, PC2)) +
  geom_hline(yintercept = 0, linetype = 2) +
  geom_vline(xintercept = 0, linetype = 2) +
  geom_point(aes(color = Ui)) +
  ggrepel::geom_text_repel(aes(label = species)) +
  scale_color_viridis_c("Functional\nUniqueness") +
  theme_bw() +
  theme(aspect.ratio = 1)

patchwork::wrap_plots(plot_reg_distinctiveness, plot_reg_uniqueness)
```

We can clearly see that, by definition, most functionally
distinct species are on the edges of the functional space (pineapple, litchi,
mango). If we compare both plots, we realize that functional distinctiveness
doesn't always come with functional uniqueness. One strange observation is that
passion-fruit has a very functional uniqueness even though it seems quite 
in the middle of the points. This is a consequence of using the raw trait
dissimilarity matrix instead of the projected axes. We're are here visualizing
the position of passion fruit in a 2D space while it should be considered in a
5D space, so it may seem close to other points in 2D while it is further away in
reality.

As was done with `mFD` to correlate the functional axes with species' traits we
can correlate functional distinctiveness to specific traits in order to see
which traits are mainly driving distinctiveness.

Regarding local level functional originality indices, the visualization can be
more difficult to grasp and depends highly on the question. Would you rather
focus on visualizing the functional distinctiveness of one species across
communities? Compare the distribution of functional distinctiveness values
across communities?

One idea to keep in mind is that averaging functional distinctiveness per
community is exactly equal to computing functional dispersion. Functional
originality is computed on a species basis, so we should be aware that if we are
rather interested by community properties than we can compute functional 
diversity metrics which are much more appropriate.

If we take again our example of the apple we can see to what extent its
functional distinctiveness varies across baskets compared to lemon.

```{r funrar-plot-local-di}
local_di_ap <- as.data.frame(sp_local_di[, c(1, 11)])
local_di_ap$basket <- rownames(local_di_ap)
```

You can of course then compute relationship between environmental covariates and
functional originality values. Because we don't have any covariable in our fruit
baskets we will stop the plotting session here.

### 3.4.2 Plotting rarity

We can plot also the rarity indices we computed and compared them to the 
functional rarity indices.

```{r funrar-plot-reg, options}
sp_di_ri <- merge(sp_di, ri, by = "species")
sp_di_ri_ui <- merge(sp_di_ri, sp_ui, by = "species")

plot_dist_ri_reg <- ggplot(sp_di_ri, aes(distinctiveness, Ri)) +
  geom_point() +
  ggrepel::geom_text_repel(aes(label = species)) +
  labs(x = "Functional Distinctiveness", y = "Geographical Restrictedness") +
  theme_bw() +
  theme(aspect.ratio = 1)

plot_dist_ri_reg
```
On this visualization we can clearly see that pineapple is overall the most
distinct species while being quite restricted in terms of baskets. On the other
hand, apple is the most functionally common and regionally widespread fruit
(low restrictedness).

We can produce similar plot with functional uniqueness:

```{r funrar-plot-reg-ui, options}
plot_ui_ri_reg <- ggplot(sp_di_ri_ui, aes(Ui, Ri)) +
  geom_point() +
  ggrepel::geom_text_repel(aes(label = species)) +
  labs(x = "Functional Uniqueness", y = "Geographical Restrictedness") +
  theme_bw() +
  theme(aspect.ratio = 1)

plot_ui_ri_reg
```
We observe a rather similar plot for functional uniqueness than for
functional distinctiveness.

We can also try to plot local scale measurements. We first convert the full
matrix of indices into easier to work with data frames through
the `matrix_to_stack()` function.

```{r funrar-data-local}
sp_local_di_df <- matrix_to_stack(
  sp_local_di, value_col = "local_di", row_to_col = "basket",
  col_to_col = "species"
)
sp_local_si_df <- matrix_to_stack(
  si, value_col = "local_si", row_to_col = "basket", col_to_col = "species"
)

sp_local_di_si <- merge(
  sp_local_di_df, sp_local_si_df, by = c("basket", "species")
)

head(sp_local_di_si)
```

Now that we have properly formatted data we can make our local plot.

```{r funrar-plot-local}
plot_local_di_si <- ggplot(sp_local_di_si, aes(local_di, local_si)) +
  geom_point(alpha = 1/3) +
  labs(x = "Functional Distinctiveness", y = "Scarcity") +
  theme_bw() +
  theme(aspect.ratio = 1)

plot_local_di_si
```
There seems to be no correlation between functional distinctiveness
and scarcity. So locally rare fruits are not necessarily distinct.

That plot concludes our plotting section on functional rarity indices.

## Corrections {.appendix}


If you see mistakes or want to suggest changes, please 
[Create an issue](https://github.com/frbcesab/workshop-free/issues)
on the source repository.



## Reuse {.appendix}



The material of this website is licensed under Creative Commons Attribution 
[CC BY 4.0](https://creativecommons.org/licenses/by/4.0/).
Source code is available at https://github.com/frbcesab/workshop-free/.


## Citation {.appendix}



Casajus N, Grenié M, Magneville C & Villéger S (2022) 
Workshop FRB-CESAB & FREE Working Group: Functional Rarity and Diversity in Ecology.