0
votes

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
1
make an example (maybe a little matrix/array) with desired outcome.Andre Elrico
How does changing NA to 0 lead to an underestimation in cumsum?Val
Maybe the function na.locf() (zoo package) is what you're looking for. This fills NA's with the last observed value after which you can use cumsum as normal.Niek
@Val: because my values differ from 0, so the cumulative sum is smaller if all NA are changed to 0Richard
@Niek na.locf() does not work with rasters or am I using it wrong?Richard

1 Answers

0
votes

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