tue_species_ranges_and_richness.Rmd

---
title: "Species ranges and richness"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(fig.width=12, fig.height=8, eval=FALSE,
                      echo=TRUE, warning=FALSE, message=FALSE,
                      tidy = TRUE, collapse = TRUE,
                      results = 'hold')
```

# Background
Patterns of species richness and occurrence density are relevant information for study design, to communicate study results or to check data quality in any multi-species occurrence dataset. Species distributions can be estimated with a large array of methods, but in some cases a rough and quick estimation of species ranges is more suitable. The extent of occurrence based on convex hulls is used in conservation and macroecology to approximate species' ranges.and range weighted endemism

# Objectives
After this exercise you will be able to:
1. Visualize species richness patterns and occurrence record density in a quick and efficient way.
2. Estimate species distribution ranges based, and create polygons of species distributions.

# Exercises
1. Load the example data, and visualize the coordinates on a map.
2. Visualize the number of occurrence records and the number of species in a lat/lon grid.
3. Estimate the range for each species using a convex hull. Write the range maps to the working directory as ESRI-shape files.
4. Visualize species richness per grid cell based on the range polygons.
5. Visualize species richness in an equal area raster.
6. Write the result rasters to a file on your disk


# Possible questions for your project
* What is the maximum range size for a species in your group
* What does the distribution of range sizes look like? Is it normally distributed?

# Tutorial
```{r}
library(tidyverse)
library(speciesgeocodeR)
library(raster)
library(rgdal)
```

## 1. Load the example data, and visualize the coordinates on a map.

```{r}
dat <- read_csv("example_data/day2_data_quality_and_sampling_bias/cleaned_occurrences_bombacoideae.csv")

# Visualize
world.inp  <- map_data("world")

ggplot()+
  geom_map(data=world.inp, map=world.inp, aes(x=long, y=lat, map_id=region), fill = "grey80")+
  xlim(min(dat$decimalLongitude, na.rm = T), max(dat$decimalLongitude, na.rm = T))+
  ylim(min(dat$decimalLatitude, na.rm = T), max(dat$decimalLatitude, na.rm = T))+
  geom_point(data = dat, aes(x = decimalLongitude, y = decimalLatitude),
             colour = "darkblue", size = 1)+
  coord_fixed()+
  theme_bw()+
  theme(axis.title = element_blank())
```

## 2. Visualize the number of occurrence records and the number of species in a grid (`RichnessGrid`)
For a quick overview you can first visualize pattern in a lat/long grid.

```{r}
# The number of occurrences per grid cell
gocc <- RichnessGrid(dat, type = "abu")

plot(gocc)

# Or more sofisticated plotting
# Prepare for plotting
plo_gocc <- data.frame(rasterToPoints(gocc))

# Visualize
ggplot()+
  geom_map(data = world.inp, map = world.inp, aes(x = long, y = lat, map_id = region), fill = "grey80")+
  xlim(min(dat$decimalLongitude, na.rm = T), max(dat$decimalLongitude, na.rm = T))+
  ylim(min(dat$decimalLatitude, na.rm = T), max(dat$decimalLatitude, na.rm = T))+
  geom_raster(data = plo_gocc, aes(x = x, y = y, fill = layer))+
  #scale_fill_viridis(name = "Species", direction = -1)+
  coord_fixed()+
  theme_bw()+
  ggtitle("Number of occurrences")+
  theme(axis.title = element_blank())

# The number of species per grid cell
gri <- RichnessGrid(dat)

plot(gri)
```

## 3. Convex hull ranges per species
Species occurrence records can be sparse for many groups, which is a problem to estimate species richness and for the bioregionalization. As the simplest model of species distribution we can use convex hull to fill the gaps. This will only work for species which range is smaller than 180 degrees longitude.

```{r}
#restrict to SOuth American and African occurrences which are the native range, except for Adansonia gregorii which is native to Autralia
add <- dat %>% 
  filter(species == "Adansonia gregorii")

dat_native <- dat %>% 
  filter(decimalLongitude > -110 & decimalLongitude < 60) %>% 
  filter(species != "Adansonia gregorii")
  bind_rows(add)

# Calculate ranges
ranges <- CalcRange(dat_native)

# Visualize
plot(ranges)

# Write to working directory
writeOGR(ranges, dsn = "example_data/day3_historical_biogeography", layer = "range_polygons_bombacoideae", 
         driver = "ESRI Shapefile", overwrite_layer = TRUE)

```

## 4. Visualize species richness per grid cell based on the range polygons.

```{r}
# Create richness raster
r_ri <- RangeRichness(ranges, res = 0.5, terrestrial = TRUE)

plot(r_ri)
```

## 5. Visualise species richness in an equal area raster
So far we have displayed species richness in a lat/lon projected raster. This is OK to get a feeling for the data, but strictly speaking wrong, since the actuals area of a 1x1 degree cell varies with latitude. So, it is better to use an equal area projection, such as the Lamberts projection.
```{r}
# Define projections
wgs1984 <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
behr <- CRS('+proj=cea +lon_0=0 +lat_ts=30 +x_0=0 +y_0=0 +datum=WGS84 +ellps=WGS84 +units=m +no_defs')

# reproject the occurrence records
pts <- dat[, c("decimalLongitude", "decimalLatitude")]%>%
  SpatialPoints(proj4string = wgs1984)%>%
  spTransform(behr)

# create a aqual area template in behrman projection
be <- raster(ncol = 360, nrow = 142, 
             xmn = -17367529, xmx = 17367529, 
             ymn = -6356742, ymx = 7348382,
             crs = behr)
be <- crop(be, extent(pts))

pts <- data.frame(dat$species, coordinates(pts))

# Equal area Occurrence number
eq_occgri <- RichnessGrid(x= pts, ras = be, type = "abu")
plot(eq_occgri)

# Equal area speces richness
eq_rigri <- RichnessGrid(x = pts, ras = be)
plot(eq_rigri)

# Range based species richness
ranges <- spTransform(ranges, behr)

eq_rri <- RangeRichness(x = ranges, ras = be)

plot(eq_rri)
```

# Output generated

- A map of the number of occurrences in a lat/lon grif
- A map of species richness in lat/lon
- A shape file with species ranges
- Plots of richness in an equal area raster.