rrastergeography

How could I add bodies of water to raster data?

I have a raster layer, which has latitudes and longitudes at a given resolution that determines the total number of grid cells. My goal is to find the distance from each cell to a body of water. How can I go about incorporating bodies of water into the data? Or, should I simply take the raster and add it to something else that already has that data? Is there a source for this? Thanks!

By the way, I did try Googling this, but all of the questions I saw are much more complicated. I've shared code below.

EDIT:

dat<-load('path/data.Rdata')
dat
class      : RasterBrick 
dimensions : 229, 216, 49464, 52  (nrow, ncol, ncell, nlayers)
resolution : 0.508, 0.24  (x, y)
extent     : -149.975, -40.247, 19.92736, 74.88736  (xmin, xmax, ymin, 
ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 
source     : memory
names      : mallar3_p//.1.1.1.1.1, mallar3_p//.2.1.1.1.1, 
mallar3_p//.1.2.1.1.1, mallar3_p//.2.2.1.1.1, mallar3_p//.1.1.2.1.1, 
mallar3_p//.2.1.2.1.1, mallar3_p//.1.2.2.1.1, mallar3_p//.2.2.2.1.1, 
mallar3_p//.1.1.1.2.1, mallar3_p//.2.1.1.2.1, mallar3_p//.1.2.1.2.1, 
mallar3_p//.2.2.1.2.1, mallar3_p//.1.1.2.2.1, mallar3_p//.2.1.2.2.1, 
mallar3_p//.1.2.2.2.1, ... 
min values :                     0,                     0,                     
0,                     0,                     0,                     0,                     
0,                     0,                     0,                     0,                     
0,                     0,                     0,                     0,                     
0, ... 
max values :          0.0007933783,          0.0010485946,          
0.0006397116,          0.0006564575,          0.0006783239,          
0.0007855534,          0.0008396646,          0.0007511564,          
0.0005148308,          0.0006255456,          0.0006670767,          
0.0006578597,          0.0006157646,          0.0003092461,          
0.0005634251, ...
Solution

  • Given your RasterBrick has 52 layers, and your intended analysis is not clear, here is a generalisable approach that will create a raster object containing distance in metres to the nearest water body. Methodologically speaking, it makes sense in your case to have a separate raster with distance values which you can then use as the basis of your analysis.

    This reprex uses lakes data from the rnaturalearth package. rnaturalearth also has river data etc. You can make a composite sf with any number of source data if you need.

    As the CRS of your RasterBrick is WGS84/EPSG:4326, you need to ensure your other data have the same CRS. In this case, WGS84/EPSG:4326 is the default CRS for rnaturalearth data.

    The cells in the distance to water body RasterLayer created by the following code align with your original RasterBrick object. This may or may not be the best approach for your analysis, but it does keep things simple. Also, as the closest water body to a given cell may be outside the extent of your RasterBrick, the extend() function is used to pad the water body RasterLayer by n number of cells (10 in this case). You can increase the extend() values if necessary.

    One further point is that raster() and sp() have been deprecated, so to future-proof your work, I encourage you to migrate to the terra package.

    library(rnaturalearth)
library(sf)
library(dplyr)
library(raster)
library(ggplot2)

# Download North America lakes and simplify
lakes_na <- ne_download(scale = 10,
                        type = "lakes_north_america",
                        category = "physical") |> 
  dplyr::select() |>
  st_combine() |>
  st_as_sf()

# Define the metadata variables from your example
nrow <- 229
ncol <- 216
nlayers <- 52
xmin <- -149.975
xmax <- -40.247
ymin <- 19.92736
ymax <- 74.88736

# Create a template for a water bodies RasterLayer based on your metadata
r <- raster(nrow = nrow, ncol = ncol, 
            xmn = xmin, xmx = xmax, ymn = ymin, ymx = ymax, 
            crs = "+proj=longlat +datum=WGS84 +no_defs")

# Extend template RasterLayer to account for water bodies near extent of your data
r <- extend(r, c(10, 10))

# Rasterize water bodies (lakes_na)
rw <- raster::rasterize(lakes_na, r, field = 1, background = NA)

# Calculate distance to the nearest water body
d2w <- distance(rw)

# Convert to df (for illustrative purposes)
d2w_df <- as.data.frame(d2w, xy = TRUE)

# Create polygon sf of extent of your RasterBrick (for illustrative purposes)
sf_poly <- st_as_sfc(st_bbox(c(xmin = xmin,
                               ymin = ymin,
                               xmax = xmax,
                               ymax = ymax))) |>
  st_as_sf(crs = 4326)

ggplot() + 
  geom_tile(data = d2w_df,
            aes(x = x, y = y, fill = layer)) +
  geom_sf(data = sf_poly,
          aes(colour = "RasterBrick extent"),
          fill = NA,
          linewidth = 1) +
  scale_fill_continuous(name = "Distance to water (m)") +
  scale_colour_discrete(name = "")

    Note that d2w extends beyond the extent of your RasterBrick: 1