I have around 80 raster stacks with 24 layers each. I'm trying to calculate the mean of all layers for each stack (e.g. mean pixel-wise in 24 layers).
I found two approaches to accomplish this but I get slightly different results in the min-max values. Here is an example:
library(raster)
set.seed(9999)
##create example raster list
rasterls=list()
for (i in 1:24){
r <- raster(nrow=(255*5), ncol=(487*5), vals=runif(255*5*487*5))
rasterls[[i]] <- r
}
##create raster stack
stackRas <- stack(rasterls)
##calculate mean using mean()
mean(stackRas, na.omit=T)
'''class : RasterLayer
dimensions : 1275, 2435, 3104625 (nrow, ncol, ncell)
resolution : 0.1478439, 0.1411765 (x, y)
extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs
source : memory
names : layer
values : 0.2453108, 0.8256463 (min, max)
'''
##calculate mean using calc(fun=mean)
calc(stackRas, fun = mean, na.omit=T)
'''
class : RasterLayer
dimensions : 1275, 2435, 3104625 (nrow, ncol, ncell)
resolution : 0.1478439, 0.1411765 (x, y)
extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs
source : C:/Users/--/AppData/Local/Temp/RtmpoVVABX/raster/r_tmp_2021-04-28_120400_4844_87419.grd
names : layer
values : 0.2138654, 0.8183816 (min, max)
'''
What is the difference between both approaches? Why the different min-max? Am I missing something? Which one should I use? I'd appreciate any insight. Thanks
There is no na.omit
argument for mean
. The correct argument should be na.rm
. When you passed a value for na.omit, it was ignored and the two methods differ in their default handling of missing values. If we use the correct argument, the results are the same.
mean(stackRas, na.rm=T)
# class : RasterLayer
# dimensions : 1275, 2435, 3104625 (nrow, ncol, ncell)
# resolution : 0.1478439, 0.1411765 (x, y)
# extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
# crs : +proj=longlat +datum=WGS84 +no_defs
# source : memory
# names : layer
# values : 0.2113316, 0.7704163 (min, max)
calc(stackRas, mean, na.rm=T)
# class : RasterLayer
# dimensions : 1275, 2435, 3104625 (nrow, ncol, ncell)
# resolution : 0.1478439, 0.1411765 (x, y)
# extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
# crs : +proj=longlat +datum=WGS84 +no_defs
# source : /tmp/RtmpRMo1uo/raster/r_tmp_2021-04-28_124738_12579_61131.grd
# names : layer
# values : 0.2113316, 0.7704163 (min, max)