[GIS] Mean by month on R stacked raster

aggregationrrasterstack

I work with MODIS NDVI rasters in 2016. I have 23 rasters stacked in one object. I have 2 raster by month. I would like the average by months and conserve a raster for each month.

ndvi.stack <- stack(result)

# attribute name for each raster stacked
idx <- seq(as.Date('2016-01-17'), as.Date('2017-01-03'), '16 day')
names(ndvi.stack) <- idx

dim(ndvi.stack)
#[1] 302 268  23
## Set up color gradient with 100 values between 0.0 and 1.0
breaks <- seq(0, 1, by=0.01)
cols <- colorRampPalette(c("red", "yellow", "lightgreen"))(length(breaks)-1)

##plot
levelplot(ndvi.stack,at=breaks, col.regions=cols, main="NDVI 2016")

I would like to get something like that

Best Answer

I have found on stack overflow a more generic way with the raster package using stackApply().

#get the date from the names of the layers and extract the month
indices <- format(as.Date(names(ndvi.stack), format = "X%Y.%m.%d"), format = "%m")
indices <- as.numeric(indices)

#sum layers
MonthNDVI<- stackApply(ndvi.stack, indices, fun = mean)
names(MonthNDVI) <- month.abb

## Set up color gradient with 100 values between 0.0 and 1.0
breaks <- seq(0, 1, by=0.01)
cols <- colorRampPalette(c("red", "yellow", "lightgreen"))(length(breaks)-1)
levelplot(MonthNDVI,at=breaks, col.regions=cols)

Et voilà

Related Solutions

[GIS] Max cell value from a stacked raster

The following example shows two ways to get at the max raster value in a stack. The first utilizes max() which also gives you a host of other useful information. The second method uses maxValue(), which gives just the max value of both the rasters in the stack

library(raster)  

# Generate some georeferenced raster data
x = matrix(rnorm(400),20,20)
rast = raster(x)
extent(rast) = c(36,37,-3,-2)
projection(rast) = CRS("+proj=longlat +datum=WGS84")

y = matrix(rnorm(400),20,20)
rast2 = raster(y)
extent(rast2) = c(36,37,-3,-2)
projection(rast2) = CRS("+proj=longlat +datum=WGS84")

raster = stack(rast, rast2)

# Now run the statistics
max(raster) # Provides min, max and additional details  # Example 1

maxValue(raster)  # Gives both values                   # Example 2...
maxValue(raster)[[1]] # Gives first in stack max value
maxValue(raster)[[2]] # Gives second in stack max value

> maxValue(raster)  # Gives both values
[1] 2.688376 2.971443
> maxValue(raster)[[1]] # Gives first in stack max value
[1] 2.688376
> maxValue(raster)[[2]] # Gives second in stack max value
[1] 2.971443
> 
> max(raster) # Provides min, max and additional details
class       : RasterLayer 
dimensions  : 20, 20, 400  (nrow, ncol, ncell)
resolution  : 0.05, 0.05  (x, y)
extent      : 36, 37, -3, -2  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
data source : in memory
names       : layer 
values      : -1.457908, 2.971443  (min, max)

[GIS] Importing multiple stacked raster images in R

Have a look at nlayers(s). The returned number of layers will equal 28 - at least for the above example with 4 multi-layer objects encompassing 7 layers each. Applying stack to multiple multi-layer files results in one huge 'RasterStack' object, i.e. all the single multi-layer objects are appended to one another.

If you would like to have separate stacks for each file, I would recommend using

s <- lapply(raster_data, stack)

which results in a list of 'RasterStack' objects, each including 7 layers rather than one huge stack. You may then access particular layers, e.g. the 2nd layer of the 3rd 'RasterStack' object, by

s[[3]][[2]]

Best Answer

Related Solutions

[GIS] Max cell value from a stacked raster

[GIS] Importing multiple stacked raster images in R

Related Question