How can I plot and overlay density plot of several rasters?

Question

I have several rasters (effectively resulted calculated from digital elevation models) (created using the raster package) and I would like to compare the distribution of the values. One possible way is to plot a density plot and overlay them on a chat.

I realise that making a stack and then the lattice package, but the rasters have different resolutions (the purpose is to test the effect of different resolutions of the resulting calculation.

I have recently become a little comfortable with the ggplot2 package, but it doesn't as I understand it, work with raster data types.

Can anyone please advise a package or technique to plot multiple density plots (and even other types such as box and whisker plots) for comparing the characteristics of different rasters please?

Answer 1

I'm not quite sure what you are after, but I'm going to take a shot anyway.

Suppose we have the following list of rasters of different resolutions and we're interested in plotting the distributions of the values inside the raster with the ggplot2 package.

library(raster)
#> Loading required package: sp
library(ggplot2)

rasterlist <- list(
  raster1 = raster(matrix(runif(100), 10)),
  raster2 = raster(matrix(rnorm(25), 5)),
  raster3 = raster(matrix(rpois(64, 2), 8))
)

What we would have to do is get this raster data into a format that ggplot2 understand, which is the long format (as opposed to wide data). This means that every observation, in the raster case every cell, should be on their own row in a data.frame. We do this by transforming each raster with as.vector() and indicate the raster of origin.

df <- lapply(names(rasterlist), function(i) {
  data.frame(
    rastername = i,
    value = as.vector(rasterlist[[i]])
  )
})
df <- do.call(rbind, df)

Now that the data is in the correct format, you can feed it to ggplot. For density plots, the x positions should be the values. Setting the fill = rastername will automatically determine grouping.

ggplot(df, aes(x = value, fill = rastername)) +
  geom_density(alpha = 0.3)

For box- or violin-plots, the group is typically on the x-axis and the values are mapped to the y-axis.

ggplot(df, aes(x = rastername, y = value)) +
  geom_boxplot()

ggplot(df, aes(x = rastername, y = value)) +
  geom_violin()

^{Created on 2020-10-04 by the reprex package (v0.3.0)}

Hope that this is roughly what you were looking for.