--- title: "belg: Boltzmann Entropy of a Landscape Gradient" author: "Jakub Nowosad" date: "`r Sys.Date()`" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{belg: Boltzmann Entropy of a Landscape Gradient} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r setup, include = FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) library(ggplot2) library(rasterVis) theme_clean = function (base_size = 12, base_family = ""){ theme_grey(base_size = base_size, base_family = base_family) %+replace% theme(axis.line = element_blank(), axis.text.x = element_blank(), axis.text.y = element_blank(), axis.ticks = element_blank(), axis.title.x = element_blank(), axis.title.y = element_blank(), axis.ticks.length = grid::unit(0, "lines"), legend.position = "none", panel.background = element_blank(), panel.border = element_blank(), panel.grid.major = element_blank(), panel.grid.minor = element_blank(), panel.spacing = grid::unit(0, "lines"), plot.background = element_blank(), strip.background = element_rect(colour = "#ffffff", fill = "#eeeeee")) } theme_set(theme_clean()) ``` Boltzmann entropy (also called configurational entropy) has been recently adopted to analyze entropy of landscape gradients (Gao et al. (2017), Gao et al. (2018)). The goal of **belg** is to provide an efficient C++ implementation of this method in R. It also extend the original idea by allowing calculations on data with missing values. # Basic example ```{r, message=FALSE} library(raster) library(belg) complex_land = raster(system.file("raster/complex_land.tif", package = "belg")) simple_land = raster(system.file("raster/simple_land.tif", package = "belg")) ``` Let's take two small rasters - `complex_land` representing a complex landscape and `simple_land` representing a simple landscape. ```{r, fig.height=4, echo=FALSE} gplot(complex_land) + geom_tile(aes(fill = value), color = "black", size = 2) + geom_text(aes(label = value), size = 4, color = "black") + coord_equal() + scale_fill_distiller(palette = "RdYlBu") + labs(title = "Complex landscape") gplot(simple_land) + geom_tile(aes(fill = value), color = "black", size = 2) + geom_text(aes(label = value), size = 4, color = "black") + coord_equal() + scale_fill_distiller(palette = "RdYlBu") + labs(title = "Simple landscape") ``` The main function in this package, `get_boltzmann()`, calculates the Boltzmann entropy of a landscape gradient: ```{r} get_boltzmann(complex_land, method = "hierarchy") get_boltzmann(simple_land, method = "hierarchy") ``` The results, unsurprisingly, showed that the complex landscape has a larger value of the Boltzmann entropy than the simple one. The `get_boltzmann()` function accepts a `RasterLayer`, `RasterStack`, `RasterBrick`, `matrix`, or `array` object as an input. As a default, it uses a logarithm of base 10 (`log10`), however `log` and `log2` are also available options for the `base` argument. ```{r} get_boltzmann(complex_land, method = "hierarchy") # log10 get_boltzmann(complex_land, method = "hierarchy", base = "log") get_boltzmann(complex_land, method = "hierarchy", base = "log2") ``` It also allows for calculation of the relative (the `relative` argument equal to `TRUE`) and absolute Boltzmann entropy of a landscape gradient. # Relative Boltzmann entropy of a landscape gradient The main idea behind the Boltzmann entropy of a landscape gradient is to calculate an entropy in a sliding window of 2 x 2 pixels. The relative configurational entropy is a sum of entropies for all windows of the original data. ```{r} get_boltzmann(complex_land, method = "hierarchy", relative = TRUE) ``` # Absolute Boltzmann entropy of a landscape gradient It is possible to calculate an average value for each sliding window of 2 x 2 pixels and therefore create a resampled version of the original dataset: ```{r, fig.height=4, echo=FALSE} complex_land_l1 = complex_land raster_template = raster(ncols = 7, nrows = 5, xmn = 0, xmx = 7, ymn = 0, ymx = 5) complex_land_l2 = raster(matrix(c(53, 32, 50, 53, 32, 81, 58, 69, 81, 58, 76, 80, 79, 76, 80, 44, 59, 62, 44, 59, 41, 47, 56, 41, 47, 57, 50, 69, 57, 50, 51, 49, 59, 51, 49), ncol = 7), template = raster_template) gplot(complex_land_l1) + geom_tile(aes(fill = value), color = "black", size = 2) + geom_text(aes(label = value), size = 4, color = "black") + coord_equal() + scale_fill_distiller(palette = "RdYlBu", limits = c(12, 98)) + labs(title = "Original dataset (Level 1)") gplot(complex_land_l2) + geom_tile(aes(fill = value), color = "black", size = 2) + geom_text(aes(label = value), size = 4, color = "black") + coord_equal() + scale_fill_distiller(palette = "RdYlBu", limits = c(12, 98)) + labs(title = "Resampled dataset (Level 2)") ``` The absolute configurational entropy is a sum of relative configurational entropies for all levels, starting from the original data to the resampled dataset with at least two rows or columns. ```{r} get_boltzmann(complex_land, method = "hierarchy", relative = FALSE) ``` # Calculation of the configurational entropy in a sliding window Determining the number of microstates belonging to a defined macrostate in a crucial concept for calculation of the configurational entropy. We explore this topic using five different cases of 2 x 2 windows: ```{r} win_1 = raster(matrix(c(1, 3, 3, 4), ncol = 2)) win_2 = raster(matrix(c(1, 3, 3, NA), ncol = 2)) win_3 = raster(matrix(c(1, 3, NA, NA), ncol = 2)) win_4 = raster(matrix(c(1, NA, NA, NA), ncol = 2)) win_5 = raster(matrix(c(NA, NA, NA, NA), ncol = 2)) ``` ## Data without missing values The configurational entropy for data without missing values is calculated using the analytical method by Gao et al. (2018). ```{r, echo = FALSE} val_cols = c("1" = "#fc8d59", "3" = "#ffffbf", "4" = "#91bfdb") gplot(win_1) + geom_tile(aes(fill = as.factor(value)), color = "black", size = 2) + geom_text(aes(label = value), size = 6) + scale_fill_manual(values = val_cols) + labs(title = "Window 1") ``` Twenty-four different microstate are possible in the above case. The common (base 10) logarithm of 24 is equal to `r round(log10(24), 6)`. We can compare this result to the `get_boltzmann()` output: ```{r} get_boltzmann(win_1, method = "hierarchy") ``` The generalized (resampled) version of this window has one value, `3`, which is a rounded average of the four original values. ## Data with missing values The papers of Gao et al. (2017, 2018) only considered data without missing values. However, the **belg** package provides a modification allowing for calculation also for data with missing values. Cells with `NA` are not considered when calculating microstates. ```{r, echo = FALSE} gplot(win_2) + geom_tile(aes(fill = as.factor(value)), color = "black", size = 2) + geom_text(aes(label = c(1, 3, 3, "NA")), size = 6) + scale_fill_manual(values = val_cols) + labs(title = "Window 2") ``` For example, three microstates are possible for the above case: ```{r, echo = FALSE} win_2_1 = win_2 win_2_2 = raster(matrix(c(3, 1, 3, NA), ncol = 2)) win_2_3 = raster(matrix(c(3, 3, 1, NA), ncol = 2)) gplot(win_2_1) + geom_tile(aes(fill = as.factor(value)), color = "black", size = 2) + geom_text(aes(label = c(1, 3, 3, "NA")), size = 6) + scale_fill_manual(values = val_cols) + labs(title = "Window 2", subtitle = "Microstate I") gplot(win_2_2) + geom_tile(aes(fill = as.factor(value)), color = "black", size = 2) + geom_text(aes(label = c(3, 3, 1, "NA")), size = 6) + scale_fill_manual(values = val_cols) + labs(title = "Window 2", subtitle = "Microstate II") gplot(win_2_3) + geom_tile(aes(fill = as.factor(value)), color = "black", size = 2) + geom_text(aes(label = c(3, 1, 3, "NA")), size = 6) + scale_fill_manual(values = val_cols) + labs(title = "Window 2", subtitle = "Microstate III") ``` The common (base 10) logarithm of 3 is equal to `r round(log10(3), 6)`. ```{r} get_boltzmann(win_2, method = "hierarchy") ``` The generalized (resampled) version of this window is 2. ```{r, echo = FALSE} win_3_1 = win_3 win_3_2 = raster(matrix(c(3, 1, NA, NA), ncol = 2)) gplot(win_3_1) + geom_tile(aes(fill = as.factor(value)), color = "black", size = 2) + geom_text(aes(label = c(1, "NA", 3, "NA")), size = 6) + scale_fill_manual(values = c("#fc8d59", "#ffffbf", "#91bfdb")) + labs(title = "Window 3", subtitle = "Microstate I") gplot(win_3_2) + geom_tile(aes(fill = as.factor(value)), color = "black", size = 2) + geom_text(aes(label = c(3, "NA", 1, "NA")), size = 6) + scale_fill_manual(values = c("#fc8d59", "#ffffbf", "#91bfdb")) + labs(title = "Window 3", subtitle = "Microstate II") ``` The third window has two combinations. The common logarithm of 2 is equal to `r round(log10(2), 6)`. ```{r} get_boltzmann(win_3, method = "hierarchy") ``` The generalized (resampled) version of this window is also 2. ```{r, echo = FALSE} gplot(win_4) + geom_tile(aes(fill = as.factor(value)), color = "black", size = 2) + geom_text(aes(label = c(1, "NA", "NA", "NA")), size = 6) + scale_fill_manual(values = c("#fc8d59", "#ffffbf", "#91bfdb")) + labs(title = "Window 4") ``` The fourth window has only one microstate, therefore its common logarithm equals to `r round(log10(1), 6)`. ```{r} get_boltzmann(win_4, method = "hierarchy") ``` The generalized (resampled) version of this window is the same as only existing value - 1. ```{r, echo = FALSE} gplot(win_5) + geom_tile(fill = "white", color = "black", size = 2) + geom_text(aes(label = c("NA", "NA", "NA", "NA")), size = 6) + labs(title = "Window 5") ``` Finally, the last window consists of four missing values. In these cases, the configurational entropy is zero. ```{r} get_boltzmann(win_5, method = "hierarchy") ``` Importantly, the generalized version of this window is represented by NA. ## References - Gao, Peichao, Hong Zhang, and Zhilin Li. "An efficient analytical method for computing the Boltzmann entropy of a landscape gradient." Transactions in GIS (2018). - Gao, Peichao, Hong Zhang, and Zhilin Li. "A hierarchy-based solution to calculate the configurational entropy of landscape gradients." Landscape Ecology 32(6) (2017): 1133-1146.