How to calculate Max value of each group used two raster stack in R?

Question

I want to apply calculations (sum) on a raster stack (contains multiple layers) according to groups of another raster stack (the groups identification is variable from one pixel time-series to another). Here is my example:

library(raster)
#1. Create raster stack "values" (rs)
set.seed(193473)
r1 <- raster(nrows = 10, ncols = 10)
r2=r3=r4=r5=r6=r7=r8=r9=r10=r11=r12=r13=r1

r1[] <- rbinom(ncell(r1), 60, prob = .1)
r2[] <- rbinom(ncell(r1), 50, prob = .1)
r3[] <- rbinom(ncell(r1), 70, prob = .1)
r4[] <- rbinom(ncell(r1), 90, prob = .1)

r5[] <- rbinom(ncell(r1), 70, prob = .1)
r6[] <- rbinom(ncell(r1), 60, prob = .1)
r7[] <- rbinom(ncell(r1), 80, prob = .1)
r8[] <- rbinom(ncell(r1), 90, prob = .1)

r9[] <- rbinom(ncell(r1), 60, prob = .1)
r10[] <- rbinom(ncell(r1), 70, prob = .1)
r11[] <- rbinom(ncell(r1), 90, prob = .1)
r12[] <- rbinom(ncell(r1), 70, prob = .1)
r13[] <- rbinom(ncell(r1), 50, prob = .1)
  
rs <- stack(r1, r2, r3, r4,r5,r6,r7,r8,r9,r10,r11,r12,r13)
nlayers(rs)
plot(rs[[1]])
plot(rs)

#2, rs_flag

rc1=function(x1) {
  ifelse(x1>7,1,0)
}

rs_flag=overlay(rs,fun=rc1)
nlayers(rs_flag)
plot(rs_flag)

#3,  rs_id

rc3= function(x3) {
  rep(seq(1,length(rle(x3)$lengths)), rle(x3)$lengths)
}

rs_id=overlay(rs_flag,fun=rc3)

plot(rs_id)

How to extract the max value of (rs) raster stack in each group of raster stack (rs_id) and replace the others value of the same group by this max value ?

I explain my problematic, see below an example of time series of one pixel :

#extract values of one pixel ex. (x,y)=(4,7)
rsp=rs[4,7]
rsp

#rsp_id
rsp_id=rs_id[4,7]
rsp_id

Extract max value of each group of stack raster (rs_id) and replace the others value of the same group by this max value.

#table preparation  
tab1 <- t(rbind(rsp, rsp_id))
tab1

#rsp_max : max value of each group

library(dplyr)
   
tab2=tab1 %>%
     group_by(rsp_id) %>%
     mutate(rsp_max=max(rsp))

tab2

At the end I get the result below for one pixel (see table below)

My question : How to calculate "rs_max" (raster stack) applying the function above between two raster stack (rs and rs_id), using "group_by" and "max" inside a function

#I try this function but not work!! :

##4# rs_id

library(dplyr)

rc4= function(x,y) {
  group_by(x) %>%
   max(y)
}



#Error in (function (x, fun, filename = "", recycle = TRUE, forcefun = FALSE,  : 
#cannot use this formula, probably because it is not vectorized

Thank you in advance for your help !

Robert Hijmans Robert Hijmans · Accepted Answer · 2020-10-01T00:14:04

I think this will do the trick

ff <- function(v) {
    x <- v[1:13]
    y <- v[14:26]
    a <- tapply(x, y, max) 
    b <- as.integer(names(a))
    as.vector(a[match(y, b)])
}

ff(as.vector(tab1))

rs_max <- calc(stack(rs, rs_id), fun=ff)
rs_max[4,7]

How to calculate Max value of each group used two raster stack in R?

1 Answers