0

I have a stack with large number of layers, and I want to run a pixel-wise cumulative sum function across the layers with reset when sum becomes <0. Ive compiled the following code:

# given that "A" is an input stack with lengthy number of layers

B<-A 

for(i in 2:nlyr(B)){
  B[[i]]<-ifel((A[[i]]+B[[i-1]])>0, A[i]]+B[[i-1]], 0 )
  
}

, which works well but I have a large stack of rasters (in terms of nlyrs). Im wondering if there is a faster way to implement this task.

tabumis
  • 106
  • 7

3 Answers3

1

UPDATED. Apparently the function I posted before was incorrect. Below is a solution which benefits from parallelization of the operation via cores argument in the terra::app function. Though I assume this is not the most optimal suggestion.

library(terra)
v <-  c(1, 2, -40, 35, 10, 10)
r <- rast(ncol=1, nrow=1, nlyr=6, vals=v)
r

cumsumreset <- function(x) {
  for (i in 2:length(x)) {
    if (x[i-1] < 0) x[i-1] <- 0
    x[i] <- x[i] + x[i-1]
    if (x[i] < 0) x[i] <- 0
  }
  x
}


a <- app(r, cumsumreset, cores=5)

a

P.S. My appologies and thanks to @RobertHijmans

tabumis
  • 106
  • 7
  • 1
    Post the correct code in the code block. The fact that ChatGPT produces broken code isn't the primary interesting thing for an answer to say; if you want to mention that ChatGPT produced a version of this with the wrong `terra::` function name, only do so after the actual answer. (And only as long as you've checked that this code truly does do what you want in all corner cases. Or probably better not to mention it at all, since [Temporary policy: ChatGPT is banned](https://meta.stackoverflow.com/q/421831) - some of the mods delete even manually-checked answers with any hint of ChatGPT :/) – Peter Cordes Feb 05 '23 at 02:17
  • ok, I hve edited my previous post – tabumis Feb 05 '23 at 10:52
1

Example data

library(terra)
v <-  c(1, 2, -40, 35, 10, 10)
cumsum(v)
#[1]   1   3 -37  -2   8  18

r <- rast(ncol=1, nrow=1, nlyr=6, vals=v)

If you want to set the negative values to zero, you can do that like this

x <- cumsum(r)
x <- ifel(x < 0, 0, x)
x
#class       : SpatRaster 
#dimensions  : 1, 1, 6  (nrow, ncol, nlyr)
#resolution  : 360, 180  (x, y)
#extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84 
#source(s)   : memory
#names       : lyr.1, lyr.2, lyr.3, lyr.4, lyr.5, lyr.6 
#min values  :     1,     3,     0,     0,     8,    18 
#max values  :     1,     3,     0,     0,     8,    18 

But that does not reset the cumulation at zero. If that is what you want you need something more complex. Perhaps like this:

cumsumreset <- function(x) {
    if (x[1] < 0) x[1] <- 0
    while(TRUE) {
        v <- cumsum(x)
        i <- which(v < 0)
        if (length(i) == 0) break
        i <- i[1]
        x[i] <- x[i] - v[i]
    }
    v
}

cumsumreset(v)
#[1]  1  3  0 35 45 55

a <- app(r, cumsumreset)
a
#class       : SpatRaster 
#dimensions  : 1, 1, 6  (nrow, ncol, nlyr)
#resolution  : 360, 180  (x, y)
#extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84 
#source(s)   : memory
#names       : lyr.1, lyr.2, lyr.3, lyr.4, lyr.5, lyr.6 
#min values  :     1,     3,     0,    35,    45,    55 
#max values  :     1,     3,     0,    35,    45,    55 
Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • Thanks @RobertHijmans. Your suggestion looks more appropriate for a bit different task: when cumulative sum should ignore negative values (by nulifying them). In my case the cumulative sum should reset to 0 upon coming across a negative value (whose absolute value is larger than the accumulated sum). This is more described here https://math.stackexchange.com/questions/4257510/mathematical-expression-of-cumulative-sum-reset-function . The updated code I posted above works well for this task – tabumis Feb 05 '23 at 21:10
  • Your use of the term "reset" is ambiguous. Apparently you just want to replace negative values with zero. That is not resetting the cumsum. – Robert Hijmans Feb 05 '23 at 21:17
  • I admit that my explanations might be ambiguous due to language barrier, will try to improve on that in my future requests – tabumis Feb 05 '23 at 22:17
  • Sorry @RobertHijmans, I was wrong – tabumis Feb 07 '23 at 23:27
1

Here is an alternative approach inspired by here Cumulative sum that resets when 0 is encountered ; this approach will reset the counter to zero

library(terra)

v <-  c(1,1,1,1,0,0,0,0,1,1,1,0,0,0,1,1)

r <- rast(ncol=1, nrow=1, nlyr=16, vals=v)

ff = function(x)
{
  cs = cumsum(x)
  cs - cummax((x == 0) * cs)
}

r2 <- app(r, ff)
plot(r2)
pdbentley
  • 430
  • 4
  • 9