-
Notifications
You must be signed in to change notification settings - Fork 4
/
regional_diversity.Rmd
88 lines (67 loc) · 2.06 KB
/
regional_diversity.Rmd
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
---
title: "Regional Diversity"
output: rmarkdown::html_vignette
vignette: >
%\VignetteIndexEntry{regional_diversity}
%\VignetteEncoding{UTF-8}
%\VignetteEngine{knitr::rmarkdown}
editor_options:
chunk_output_type: console
---
```{r, include = FALSE}
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>",
warning = FALSE,
message = FALSE)
```
```{r setup}
library(obisindicators)
library(dplyr)
library(sf)
```
## Get regional biological occurrences
Let's use the regional subset for the South Atlantic from the full OBIS dataset otherwise available at https://obis.org/data/access.
```{r occ}
occ <- occ_SAtlantic # occ_1M OR occ_SAtlantic
```
## Create a global h3 hexagonal grid
```{r h3construct}
hex_res <- 1 # hex_res 0 is too big to work, all others work
hex <- obisindicators::make_hex_res(hex_res)
# mapview::mapview(hex) # you can view the hex grid with h3 IDs
# === Then assign cell numbers to the occurrence data:
occ <- occ %>%
mutate(
cell = h3::geo_to_h3(
data.frame(decimalLatitude, decimalLongitude),
res = hex_res))
```
## Calculate indicators
The following function calculates the number of records, species richness, Simpson index, Shannon index, Hurlbert index (n = 50), and Hill numbers for each cell.
Perform the calculation on species level data:
```{r calc_indicators}
idx <- calc_indicators(occ)
```
Add cell geometries to the indicators table (`idx`):
```{r dgcellstogrid}
grid <- hex %>%
inner_join(
idx,
by = c("hexid" = "cell"))
# you could now visualize with:
# plot(grid["es"])
# mapview::mapview(grid["es"])
```
## Plot maps of indicators
Let's look at the resulting indicators in map form.
```{r map_gcs}
# ES(50)
gmap_indicator(grid, "es", label = "ES(50)")
# Shannon index
gmap_indicator(grid, "shannon", label = "Shannon index")
# Number of records, log10 scale, Robinson projection (default)
gmap_indicator(grid, "n", label = "# of records", trans = "log10")
# Number of records, log10 scale, Geographic projection
gmap_indicator(grid, "n", label = "# of records", trans = "log10", crs=4326)
```