Search code examples
rrasternacumsum

Cumulative sum: deal with NA pixel in raster data in R

I have a multiple band raster and want to calculate cumulative sums of each pixel throughout the bands. The bands include NA pixels at varying locations and the cumsum function in calc from the raster package writes NA in subsequent bands if it encounters a NA value in a pixel. Finally I have only NA pixels left in the last cumulative sum band. I changed NA to zero, but that leads to an underestimation of the values.

Is it possible to use the preliminary value to NA? Or maybe even an average of the preliminary and the following value?

library(raster)
raster <- stack("inputPath")
cumulative_sum <- calc(raster, cumsum)

here is an example of what I mean

input band 1
1    4    7
NA   5    8
3    6   NA

input band 2
2   NA   NA
3    6    9
4    7   10

input band 3
1    4    7
2    5    8
3    6    9

result with calc and cumsum
4    NA   NA
NA   16   25
10   19   NA


desired output (last resulting band <- band1 + band2+ band3)
4    12   21
5    16   25
10   19   19
Solution

  • Using @Niek 's suggestion, you can implement it this way:

    library(raster)
library(zoo)

# Generate example data

set.seed(42)
r <- raster(ncols=3,nrows=3)

r3 <- do.call(stack,replicate(3,setValues(r,sample(c(runif(ncell(r)),NA),ncell(r),replace = T))))

f <- function(ii){

  if(is.na(ii[1])) ii[1] <- 0

  cumsum(na.locf(ii))

}

r3c <- calc(r3,fun = f)

    As you can see, any NA values in layer 1 are set to 0, since there's no value to carry forward. You could also do this prior to calc and remove the if clause from f.

    To check we'll plot it and you'll see that the NA's are gone:

    plot(is.na(r3),main=paste('Input',c(1:3)))

    enter image description here

    plot(is.na(r3c),main=paste('Cumsum',c(1:3)))

    enter image description here