I have a raster stack of 7 rasters. Each raster has the following properties
class : RasterLayer
dimensions : 21600, 43200, 933120000 (nrow, ncol, ncell)
resolution : 0.008333333, 0.008333333 (x, y)
extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs
source : tmax.tif
names : tmax
values : 1, 100 (min, max)
I want to create a single raster that has for each grid, the maximum value across the 7 rasters. I thought I would do this following a solution proposed here:
select highest value in Raster stack and show layer name as legend
max_raster <- which.max(my_raster_stack)
df <- reshape2::melt(as.matrix(max_raster))
df$value <- factor(var_name[df$value], var_name)
However, due to sheer size of the raster stack, it's taking my hours to do this (still waiting for it to finish). Does anyone know any quicker way to do this?
You can use max to get the maximum value of each cell. Here with terra (the replacement of raster, but it works the same with raster)
Example data
library(terra)
f <- system.file("ex/logo.tif", package="terra")
r <- rast(f, f) * 1:6
r
#class : SpatRaster
#dimensions : 77, 101, 6 (nrow, ncol, nlyr)
#resolution : 1, 1 (x, y)
#extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs
#source : memory
#names : lyr1, lyr2, lyr3, lyr4, lyr5, lyr6
#min values : 0, 0, 0, 0, 0, 0
#max values : 255, 510, 765, 1020, 1275, 1530
Solution
x <- max(r)
x
#class : SpatRaster
#dimensions : 77, 101, 1 (nrow, ncol, nlyr)
#resolution : 1, 1 (x, y)
#extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs
#source : memory
#name : max
#min value : 0
#max value : 1530
To find out, for each grid cell, which layer has the maximum value (or the first one if there are ties), you can do
w <- which.max(r)
w
#class : SpatRaster
#dimensions : 77, 101, 1 (nrow, ncol, nlyr)
#resolution : 1, 1 (x, y)
#extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs
#source : memory
#name : which.max
#min value : 1
#max value : 6