Summary

Raster R package is really useful for geographic data analysis and modeling due to its capacity of reading, creating, manipulating and analyzing of spatial data with general raster data manipulation functions (Hijmas, 2020; Hijmas, 2021). With this package we can handle several research objectives related with remote sensing, spatial ecology, land use and land cover issues, agriculture, urban expansion, climate change, biological invasions, etc., that involves multitemporal information in raster data. In this case, we want to get habitat information from land cover data of the continental part of Colombia, as an external input data to improve species distributions models (SDM). The land cover product used here is the global Land cover classification gridded maps derived from satellite observations, and the output result is a habitat map as defined in the International Union for Conservation of Nature (IUCN) habitat classification scheme (Jung et al., 2020). The map had to be reclassified to set links between ESA land covers and IUCN habitats.

Introduction

The patterns that shape the spatial and temporal distribution of the species is one of the most fundamental pieces of ecological information (Beale & Lennon, 2012). Species distributions models (SDMs) explore the relationship between geographical occurrences of species and corresponding environmental variables (Naimi & Araujo, 2016) and essentially are based in spatial (and some in temporal) interpolation and extrapolation for the present or into the future. The accuracy of each model will depend on how can handle the different processes that a species have to occur within an ecological community like dispersal, environmental tolerance, species interactions, biogeographical history, etc.(Li et al., 2020).

Species distribution modeling has grown rapidly in the last years and with these the software and tools available to display this analysis (specially MaxEnt standalone applications and R programming packages like sdm, zoon R, phyr, etc), and select the most reliable variables for efficient calibrations (kuenm, ENMTools, ENMeval, MaxentVariableSelection, etc) (Golding et al., 2017; Cobos et al., 2019a; 2019b). Each model with its combination of parameters can identify uncertainty in different levels to see where future improvements can be made to diminish overprediction or subprediction on species distribution. One of the ways to identify sources of uncertainty is through acknowledging the species habitat.

According to the definition of Kearney (2006), habitat is a description of a physical place, at a particular scale of space and time, where an organism either actually or potentially lives. So we have an approximation of the habitat if we know the land covers within a particular area of interest, and this combined with sufficient ecological information of a specie could be an interesting input to evaluate a candidate model. This work flow is the first of two documents indicating the process to homologate land covers with habitats. This document will show the homologation between raster land cover data from European Space Agency and habitat classification scheme of International Union for Conservation of Nature (IUCN).

Package overview

The functions used here from the package raster are as follows:

  • raster() load the raster document with land cover data.
  • reclassify() the land cover data to make the homologation.
  • projectRaster() change the projection of the input raster into another raster of interest.
  • resample() transfers values between non matching Raster objects in terms of origin and resolution.
  • overlay() create a new raster based on two or more raster objects (Similar to raster calculator in SIG software)
  • rasterToPoints() Raster to point conversion.
  • rasterFromXYZ() Create a Raster object with x,y and z values, being x and y the spatial coordinates.
  • writeRaster() write an entire Raster object to a file.

Package work flow

The general work flow when using the raster package to homologate covers and habitats, starts with the load of the Raster object with land cover data. Then, we need to reclassify the land cover categories and suit them to the IUCN scheme classification. There are some important habitats that can not be only shaped with land cover data, so it has to include altitude data with a DEM or another source that you consider pertinent. Most likely, the altitude data have different attributes from the land cover map. If that is the case, the function projectRaster() and resample() will help to align both maps in pixel size, extent and dimensions.

The final map will be constructed with two different raster objects: DEM and land cover data, and some categories will need simple map algebra operations. For this we will use the function overlay(). At this point the final map of habitats is constructed, but can be modified their aspect to inspect the results before the final output been wrote into a new file.

Usage

Here we will provide an introduction to using the raster R package to make a homologation between land cover data and habitats, starting with ESA raster data.

Global land cover: download data

Land covers can be a suitable input to estimate the habitat of a species. Some products are constantly updated, which provides the advantage of temporality in the application of this routine. In this case, the product ** Land cover classification derived from satellite observations ** from the European Space Agency (ESA) will be used, which can be downloaded at the following link: https://cds.climate.copernicus.eu/cdsapp#!/dataset/satellite-land-cover?tab=overview. This product is updated annually and has a spatial resolution of 300 meters.

Load data

Depending on the study area, it may be necessary to work with more than one raster at the same time or a single raster may have to be masked with a shape. This will depend on the source you want to work with. These actions can be performed in RStudio using the merge and mask functions, or with the help of a GIS software. In this case, the study area corresponding to the continental part of Colombia has already been cut from the ESA coverage map, with the help of the QGIS software. For the meaning of the ESA product codes, go to the following link: https://datastore.copernicus-climate.eu/documents/satellite-land-cover/D3.3.12-v1.3_PUGS_ICDR_LC_v2.1.x_PRODUCTS_v1.3.pdf 

esa_col <- raster("ESA_Col_Cob_2012.tif", na.rm= TRUE)
plot(esa_col, xlim= c(-80,-65), ylim= c(-5,15), main= "Coberturas ESA - Colombia")

Reclassification Covers-Habitat

Among the elements that make up the habitat, the land covers offer an approximation to the physical space used by a species, so if the covers used by a species are recognized, it is possible to get closer to identifying its habitat. The following is a proposal for homologation between land cover and habitats according to the classification scheme of the International Union for the Conservation of Nature (IUCN). To know the meaning of the IUCN habitats, go to the following link: https://www.iucnredlist.org/resources/habitat-classification-scheme. Not all IUCN habitats were used, as it is not possible to assign each one a cover. 

mtx <- matrix(c(10,1401,11,1402,12,1401,30,1401,40,1406,50,105,60,105,61,105,80,101,90,101,100,1406,110,1406,120,305,122,305,130,405,150,600,153,600,160,107,170,107,180,406,190,1404,200,600,210,500,220,1700), ncol = 2, byrow = TRUE)
mtx
##       [,1] [,2]
##  [1,]   10 1401
##  [2,]   11 1402
##  [3,]   12 1401
##  [4,]   30 1401
##  [5,]   40 1406
##  [6,]   50  105
##  [7,]   60  105
##  [8,]   61  105
##  [9,]   80  101
## [10,]   90  101
## [11,]  100 1406
## [12,]  110 1406
## [13,]  120  305
## [14,]  122  305
## [15,]  130  405
## [16,]  150  600
## [17,]  153  600
## [18,]  160  107
## [19,]  170  107
## [20,]  180  406
## [21,]  190 1404
## [22,]  200  600
## [23,]  210  500
## [24,]  220 1700
esa_hab_rcl <- raster::reclassify(x= esa_col, rcl = mtx)
plot(esa_hab_rcl, main= "Homologación Coberturas-Habitats") #Resolucion espacial de 300 m

Add altitude variable

The IUCN has habitats with a very particular geographic delimitation for Colombia, but it is not possible to assign a coverage category to them. These habitats include 1.9. Forest - Subtropical/Tropical moist montane, 3.7 Shrubland - High altitude and 4.7 Grassland - High altitude. Many species use these habitats in Colombia, so they must be differentiated; For this, a DEM of 90 meters will be used that will allow the delimitation of these habitats. To do this correctly, the DEM must be reprojected and have the same extent, resolution and dimensions as the ESA map. 

dem_colombia <- raster("dem_colombia3.tif")
dem_proj <- projectRaster(dem_colombia, crs = crs(esa_hab_rcl))
dem_al <- resample(dem_proj, esa_hab_rcl, method= "bilinear")

It was proposed to establish a limit of 500 meters to differentiate lowland from montane habitats.

esa_hab <- overlay(dem_al, esa_hab_rcl, 
                       fun= function(altura, cover){
  ifelse((altura >= 500) & (cover == 105),109,
         ifelse((altura >= 500) & (cover == 305),307,
                ifelse((altura >= 500) & (cover == 405),407,cover)))
})

plot(esa_hab, main= "Habitats_ESA_300m")

Map visualization

This is a representation with colors more suited to the kinds of habitats that can be generated. However, as it requires an additional transformation of the data, it is suggested to work with the map esa_hab.

df <- data.frame(rasterToPoints(esa_hab))
raster_test <- rasterFromXYZ(df[1:3])
cut <- sort(unique(esa_hab))
mypal <- c("darkcyan", "cornflowerblue", "chartreuse4", "darkorchid", "red", "chartreuse", "chartreuse4", "aquamarine3", "blue3", "darkorange", "blanchedalmond", "pink", "black", "darkkhaki")
plot(raster_test,
     main= "Habitats_ESA_300m",
     col= mypal,
     breaks= cut)

Create raster file

This habitat map built from ESA covers can be an important input in updating and refining routines for the construction of species distribution models. For this reason, it is important to save a copy in .tif format in the working directory.

#writeRaster(esa_hab, "Habitats_ESA.tif", overwrite= T)

References

Beale, C. M. & Lennon, J. J. (2012). Incorporating uncertainty in predictive species distribution modelling. Philosophical Transactions of the royal society B, 367, 247-258.

Cobos, M. E., Peterson, A. T., Barve, N. & Osorio-Olvera, L. (2019a). kuenm: an R package for detailed development of ecological niche models using Maxent. Bioinformatics and genomics. PeerJ 7:e6281 https://doi.org/10.7717/peerj.6281

Cobos, M. E., Peterson, A. E., Osorio-Olvera, L. & Jiménez-García, D. (2019b). An exhaustive analysis of heuristic methods for variable selection in ecological niche modeling and species distribution modeling. Ecological Informatics, 53, 1-10.

Golding, N., August, T. A., Lucas, T. C. D., Gavaghan, D. J., van Loon, E. E. & McInerny, G. (2017). The zoon r package for reproducible and shareable species distribution modelling. Methods in ecology and evolution, 9, 260-268.

Hijmas, R. (2020). Raster package v.3.4-5. Available in https://www.rdocumentation.org/packages/raster/versions/3.4-5

Hijmas, R. (2021). The raster package. Available in https://rspatial.org/raster/pkg/RasterPackage.pdf

Jung, M., Dahal, P. R., Butchart, S. H. M., et al. (2020). A global map of terrestrial habitat types. Nature - Scientific data, 7:256, 1-8.

Kearney, M. (2006). Habitat, environment and niche: what are we modelling?. OIKOS, 115, 186-191.

Li, D., Dinnage, R., Nell, L. A., Helmus, M. R. & Ives, A. R. (2020). phyr: An r package for phylogenetic species‐distribution modelling in ecological communities. Methods in ecology and evolution, 11, 1455-1463.

Naimi, B. & Araújo, M. B. (2016). sdm: a reproducible and extensible R platform for species distribution modelling. Ecography, 39, 368-375.