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

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.

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:

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.

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.

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:

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.

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:

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).

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:

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.

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:

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).

get_boltzmann(win_2, method = "hierarchy")

The generalized (resampled) version of this window is 2.

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).

get_boltzmann(win_3, method = "hierarchy")

The generalized (resampled) version of this window is also 2.

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).

get_boltzmann(win_4, method = "hierarchy")

The generalized (resampled) version of this window is the same as only existing value - 1.

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.

get_boltzmann(win_5, method = "hierarchy")

Importantly, the generalized version of this window is represented by NA.

References



Nowosad/belg documentation built on Jan. 1, 2023, 3:27 p.m.