/
restrictedness.R
176 lines (155 loc) · 5.8 KB
/
restrictedness.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
# Functions to compute restrictedness on various datasets
#
# Authors: Pierre Denelle & Matthias Grenié
#
#
#' Geographical Restrictedness for stacked data.frame
#'
#' Compute the geographical restrictedness for each species present in the
#' stacked data.frame. Geographical restrictedness is an index related to the
#' extent of a species in a given dataset, it is close to 1 when the species is
#' present in only a single site of the dataset (restricted) and close to 0 when
#' the species is present at all sites. It estimates the geographical extent of
#' a species in a dataset. See [restrictedness()] for details
#' on restrictedness computation. You can either use `_stack()` or `_tidy()`
#' functions as they are aliases of one another.
#'
#' @param com_df a stacked (= tidy) data.frame of communities
#'
#' @param sp_col a character vector indicating the name of the species column
#'
#' @param com a character vector indicating the name of the community column
#'
#' @param relative a logical (default = FALSE), indicating if restrictedness
#' should be computed relative to restrictedness from a species occupying a
#' single site
#'
#' @return A stacked data.frame containing species' names and their
#' restrictedness value in the **Ri** column, similar to what
#' [uniqueness_stack()] returns.
#'
#' @seealso [restrictedness()], [uniqueness_stack()]
#'
#' @examples
#' data("aravo", package = "ade4")
#'
#' # Site-species matrix converted into data.frame
#' mat = as.matrix(aravo$spe)
#' dat = matrix_to_stack(mat, "value", "site", "species")
#' dat$site = as.character(dat$site)
#' dat$species = as.character(dat$species)
#' ri_df = restrictedness_stack(dat, "species", "site")
#' head(ri_df)
#'
#' @export
restrictedness_stack = function(com_df, sp_col, com, relative = FALSE) {
# Test to be sure of inputs
full_df_checks(com_df, sp_col, com)
# get the total number of communities
n_com = length(unique(com_df[, com]))
# Compute the sum of all species' occurrences divided by n_com
occupancy = table(com_df[, sp_col]) / n_com
# Format occupancy in data.frame
occupancy = data.frame("sp" = names(occupancy), "Ri" = as.numeric(occupancy))
colnames(occupancy)[1] = sp_col
occupancy$Ri = 1 - occupancy$Ri
# Standardize if suggested by R_i of species present in a single site
if (!relative) {
r_one = 1
} else {
r_one = 1 - 1 / n_com
}
occupancy$Ri = occupancy$Ri / r_one
return(occupancy)
}
# Restrictedness tidy alias
#' @export
#' @rdname restrictedness_stack
restrictedness_tidy = restrictedness_stack
#' Geographical Restrictedness on site-species matrix
#'
#' Computes geographical restrictedness from a site-species matrix.
#' Geographical restrictedness is an index related to the extent of a species
#' in a given dataset, it is close to 1 when the species is present in only a
#' single site of the dataset (restricted) and close to 0 when the species is
#' present at all sites. It estimates the geographical extent of a species in a
#' dataset. See `Details` section to have details on the formula used for
#' the computation. The sites-species matrix should have **sites**
#' in **rows** and **species** in **columns**, similar to \pkg{vegan} package
#' defaults.
#'
#' @param pres_matrix a site-species matrix, with species in rows and sites
#' in columns, containing presence-absence, relative abundances or
#' abundances values
#'
#' @param relative a logical (default = FALSE), indicating if restrictedness
#' should be computed relative to restrictedness from a species occupying a
#' single site
#'
#' @return A stacked data.frame containing species' names and their
#' restrictedness value in the **Ri** column, similar to what
#' [uniqueness()] returns.
#'
#' @details
#' Geographical Restrictedness aims to measure the regional extent of a species
#' in \pkg{funrar} it is computed the simplest way possible: a ratio of the
#' number of sites where a species is present over the total number of sites in
#' the dataset. We take this ratio off 1 to have a index between 0 and 1 that
#' represents how restricted a species is:
#' \deqn{
#' R_i = 1 - \frac{N_i}{N_tot},
#' }{
#' R_i = 1 - (N_i/N_tot),
#' }
#' where \eqn{R_i} is the geographical restrictedness value, \eqn{N_i} the total
#' number of sites where species \eqn{i} occur and \eqn{N_tot} the total number
#' of sites in the dataset.
#' When `relative = TRUE`, restrictedness is computed relatively to the
#' restrictedness of a species present in a single site:
#' \deqn{
#' R_i = \frac{R_i}{R_one}
#' }{
#' R_i = R_i / R_one
#' }
#' \deqn{
#' R_i = \frac{1 - \frac{K_i}{K_tot}}{1 - \frac{1}{K_tot}}
#' }{
#' R_i = (1 - K_i/K_tot)(1 - 1/K_tot)
#' }
#' \deqn{
#' R_i = \frac{K_tot - K_i}{K_tot - 1}
#' }{
#' R_i = (K_tot - K_i)(K_tot - 1)
#' }
#' Other approaches can be used to measure the geographical extent
#' (convex hulls, occupancy models, etc.) but for the sake of simplicity only
#' the counting method is implemented in \pkg{funrar}.
#'
#' @examples
#' data("aravo", package = "ade4")
#' # Site-species matrix
#' mat = as.matrix(aravo$spe)
#' ri = restrictedness(mat)
#' head(ri)
#'
#' @export
restrictedness = function(pres_matrix, relative = FALSE) {
# Check site-species matrix type
check_matrix(pres_matrix, "site-species")
# get the total number of communities
n_com = nrow(pres_matrix)
# Convert all the matrix values into 0/1
pres_matrix[which(is.na(pres_matrix))] = 0
pres_matrix[which(pres_matrix > 0)] = 1
# Compute the sum of all species' occurrences divided by n_com
occupancy = 1 - (colSums(pres_matrix, na.rm = TRUE) / n_com)
if (!relative) {
r_one = 1
} else {
r_one = 1 - 1 / n_com
}
# Format occupancy in data.frame
occupancy = data.frame("species" = names(occupancy),
"Ri" = as.numeric(occupancy) / r_one)
return(occupancy)
}