some of the code is based on the vignette of the dismo package (Hijmans & Elith 2011) and from our friend Vijay Barve (https://vijaybarve.wordpress.com)

So, we’ve already seen where do the primary biodiversity data comes from (e.g. point occurrences from inventories, collections, museums, etc.). Some of this data is already available in the internet!!

Exercise 1

Install and load the following packages

library(rgdal)
library(rgeos)
library(sp)
library(raster)
library(maptools)
library(maps)
library(rgbif)
library(dismo)
library(ggmap)
library(rvertnet)

Get some occurrence data for a particular species from GBIF, directly within R. This may take some time, given the number of occurrences for the selected species. Here, we’ll work with the fossa! NOTA: we need to have an internet connection

fossa <- gbif("Cryptoprocta","ferox")

Take a look at the object dimensions (i.e. how many lines and columns)

dim(fossa)

Check the column names, so we can look for the coordinates later

names(fossa)

Create another object keeping only the coordinates (longitude/latitude)

fossa_points <- fossa[,c("lon","lat")] 

NOTE: the variable names’ may be different (.e.g “LATITUDE”, “Latidude”, “lat”, etc.). Hence, it is always good to check the column names (as above)

Plot the points (x,y)

plot(fossa_points[,1],fossa_points[,2],col="blue",pch=19)

Add a world map, for reference. We need to load a dataobject from the maptools package

data(wrld_simpl)
plot(wrld_simpl,add=T)

Are there any “weird” points?

Do the opposite, first plot the worldmap and then the points

plot(wrld_simpl,col="light yellow")

Add the points

points(fossa_points[,1],fossa_points[,2],col="blue",cex=0.2)

NOTE: the arguments pch above and cex here are used simply to define the size and tipe of the points. You can play with this!

Depending on the species, the map/points can be too large or be outside our region of interest. We can focus on a particular region/country. Let’s do this for Madagascar

#Firts, get the "position" of Madagascar in the "wrld_simpl" object
which(wrld_simpl$NAME=="Madagascar")
#Then, use these indices (numbers) to subset the original object
mdg.nmb_map <- wrld_simpl[108,]

Now you can repeat the steps above (plotting) for the selected region (Madagascar) only.

Data checking & cleaning

Check if there are any duplicated points

fossa_dups <- duplicated(fossa_points)
### NOTE: the function "duplicated" returns the results of a logical test (e.g. TRUE or FALSE)
# How many are duplicates?
length(which(fossa_dups==TRUE))
# How many are NOT duplicates?
length(which(fossa_dups==FALSE))
# Keep only those lines that are not duplicates
fossa_dups_row <- which(fossa_dups==TRUE)
# What's the size? That is, how many points are duplicates
length(fossa_dups_row)
# Create another object withoyt the duplicate records
fossa_nodups <- fossa_points[-fossa_dups_row,]
# What are the dimensions of the new object?
dim(fossa_nodups)
# Take a look at the first rows of data
head(fossa_nodups)

Check if there are data with longitude/latitude of zero

fossa_lonzero = subset(fossa_nodups, lon==0)

fossa_latzero = subset(fossa_nodups, lat==0)

NOTE: Remember that variable names can be different (Latitude, latitude, lat, etc…)

Check if there are data with NAs

fossa_nas <- which(is.na(fossa_nodups[,1])==TRUE)
fossa_points_nonas  <- fossa_nodups[-fossa_nas,]

#how many NAs?
dim(fossa_nodups)
dim(fossa_points_nonas)

Exercise 2

Plot with other maps, using other data (e.g. download other species using another package: “rgbif”).

civet <- occ_data(scientificName = 'Fossa fossana', hasCoordinate = TRUE)

NOTA: You can choose any other species (try!)

Get a world map

world = map_data("world")

Plot the points using the map as a polygon and defining the points’ attributes (colour, size, etc)

ggplot(world, aes(long, lat)) +
geom_polygon(aes(group = group), fill = "white", 
              color = "gray40", size = .2) +
geom_jitter(data = civet$data,
aes(civet$data$decimalLongitude, civet$data$decimalLatitude), alpha=0.6, 
             size = 4, color = "red") +
labs(title = "Malagasy civet")

Is it good? What if we use another species?

Exercise 3

Now with GoogleMaps

ayeaye = occ_data(scientificName = 'Daubentonia madagascariensis', hasCoordinate = TRUE)

# Define the "extent", the map's "window"
e = extent(43, 51, -26, -12)
# Take that window from GoogleMaps
r = gmap(e)
# plot
plot(r, interpolate=TRUE, main="Map")
# join with the lat/long values
xy1=cbind(ayeaye$data$decimalLongitude,ayeaye$data$decimalLatitude)
# Add the points, using the same projection as the map (i.e. "Mercator")
points(Mercator(xy1) , col='blue', pch=20)
# Write the species name
text(4500000,-2400000,"\n\n\nDaubentonia \nmadagascariensis", adj = c(0,1), col="red", cex=0.75)

Exercise 4

Range maps from point data

This part can be used to create “simple” range maps based on geometry (e.g. minimum convex polygons, etc.), without considering environmental variables (no ENMs or SDMs).

Convex hull (minimum convex polygon)

This model draws a convex hull around all ’presence’ points.

#create a polygon around the species' records
fossa_hull <- convHull(fossa_points_nonas)

#Plot the created polygon
plot(mdg.nmb_map)
plot(fossa_hull,col="green",add=T)

#Create a raster to be used as a mask
r <- raster()
values(r) <- 0.5
ext <- c(43, 51, -26, -12)
r1 <- crop(r,ext)
#Use the raster and the polygon to create a new "map"
fossa_hullmodel <- predict(fossa_hull, r1)
plot(fossa_hullmodel, main='Convex Hull')
plot(mdg.nmb_map, border='dark grey',add=T) 
points(fossa_points_nonas, pch='+')

Circles

This model draws circles around all ’presence’ points.

fossa_circ <- circles(fossa_points_nonas, lonlat=TRUE)
fossa_circles <- predict(r1, fossa_circ, mask=TRUE)

plot(fossa_circles, main='Circles')
points(fossa_points_nonas, pch=19, cex=0.5, col="red")

#Play with the distance to create the "circles"
fossa_circ2 <- circles(fossa_points_nonas, d=100000, lonlat = TRUE)
fossa_circles2 <- predict(r1, fossa_circ2, mask=TRUE)

plot(fossa_circles2, main='Circles')
points(fossa_points_nonas, pch=19, cex=0.5, col="red")

NOTA: Save the R “workspace”. That is, everything that was done in these exercises must be saved in an “.Rdata” file. These objects may be used later on. DO NOT CLOSE R IF WITHOUT SAVING THE FILE!!