Function to sum each grid cells of raster stack using other rasters as an indicator
## input raster
s < stack(list.files("~/dailyraster", full.names=TRUE)) # daily raster stack
r_start < raster("~/stackSumSTART.asc") # this raster contain starting Julian day
r_end < raster("~/stackSumEND.asc") # this raster contain ending Julian day
noNAcells < which(!is.na(r[])) # cell numbers which contain values
## dummy raster
x < r
x[] < NA
## loop
for (i in noNAcells) {
x[i] < sum(s[[r_start[i]:r_end[i]]][i])
}
I would like to create a function like stackApply()
, but I want it to work on a cell basis.
Above is a for()
loop version and it works well, but it takes too much time.
The point is that each cell gets the range of sum()
from two raster layers, r_start
, r_end
in above script.
Now I am struggling to transform this code using apply()
family.
Is there any possibility to improve the speed with for()
loop? or please give me some tips to write this code in apply()
Any comments will help me, thank you.
1 answer

Your approach
x < s$layer.1 system.time( for (i in 1:ncell(x)) { x[i] < sum(s[[r_start[i]:r_end[i]]][i], na.rm = T) } )
user system elapsed 0.708 0.000 0.710
My proposal
You can add the rasters used as indices at the end of your stack and then use
calc
to highly speed up the process (~3050x).s2 < stack(s, r_start, r_end) sum_time < function(x) {sum(x[x[6]:x[7]], na.rm = T)} system.time( output < calc(s2, fun = sum_time) )
user system elapsed 0.016 0.000 0.015
all.equal(x, output)
[1] TRUE
Sample Data
library(raster) # Generate rasters of random values r1 < r2 < r3 < r4 < r5 < r_start < r_end < raster(ncol=10, nrow=10) r1[] < rnorm(ncell(r1), 1, 0.2) r2[] < rnorm(ncell(r2), 1, 0.2) r3[] < rnorm(ncell(r3), 1, 0.2) r4[] < rnorm(ncell(r4), 1, 0.2) r5[] < rnorm(ncell(r5), 1, 0.2) s < stack(r1,r2,r3,r4,r5) r_start[] < sample(1:2, ncell(r_start),replace = T) r_end[] < sample(3:5, ncell(r_end),replace = T)