Extract pixel values in raster stack using shp

Question

I am using the following code to extract mean pixel values over an area represented by a polygon:

mean<-lapply(S2_stack, FUN=function (S2_stack) {data.frame(mean=extract(S2_stack, polygons2, fun=mean))})

Where:

> S2_stack
class      : RasterStack 
dimensions : 1454, 1595, 2319130, 4  (nrow, ncol, ncell, nlayers)
resolution : 10, 10  (x, y)
extent     : 744490, 760440, 4773400, 4787940  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=32 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 
names      : Subset_S2_MSIL2A_20180511T100029_MCI, Subset_S2_MSIL2A_20180710T100029_MCI, Subset_S2_MSIL2A_20180906T101021_MCI, Subset_S2_MSIL2A_20181115T101251_MCI 

and polygons2:

polygons2 <- readOGR('/path/Insitu_poly.shp')

with:

> class(polygons2)
[1] "SpatialPolygonsDataFrame"
attr(,"package")
[1] "sp"

When I try to run my code I am getting the following error:

 Error in (function (classes, fdef, mtable)  : 
  unable to find an inherited method for function ?extract? for signature ?"numeric", "SpatialPolygonsDataFrame"?

Not sure why. It seems extract needs a SpatialPolygonsDataFrame and as far as I know, it should run

Any idea?

Answer 1

Assuming that the correct way to proceed in this case is simply to do

mean <- extract(S2_stack, polygons2, fun=mean)

as said in the comments, however, if you want to pass your stack in lapply you have to unstack it first.

mean<-lapply(unstack(S2_stack), FUN=function (S2_stack) 
{data.frame(mean=extract(S2_stack, 
polygons2, fun=mean))})

Answer 2

I recommend you to use the exactextractr package instead of the extract function of the raster package which is extremely slow. Here is an example of how to use it:

if (!require("raster")) install.packages("raster")
if (!require("rgdal")) install.packages("rgdal")
if (!require("exactextractr")) install.packages("exactextractr")

suppressPackageStartupMessages({
library(raster)
library(rgdal)
library(exactextractr)
})

I assume two variables loaded. img is a RasterStack with several layers and shp is a ESRI Shapefile containing polygons. This shapefile contains a field called class with the land cover. You can load this variables using: img <- raster::stack('file.tif') and shp <- rgdal::readOGR('file.shp').

# extract pixel values on a list of data.frames where
# length(ls.df) == number of polygons of the ESRI Shapefile
ls.df <- exact_extract(img, shp, include_cols="class", progress=TRUE)
# filter only full pixel (coverage_fraction==1)
for(i in 1:length(ls.df)){
ls.df[[i]] <- ls.df[[i]][ls.df[[i]]$coverage_fraction == 1,]
}
# create a single data.frame containing all the pixel values
df <- data.frame()
for(i in 1:length(ls.df)){
df <- rbind(df, ls.df[[i]])
}

The exact_extract function computes the coverave fraction of each pixel covered by the polygon. You remove this information using as follow:

# remove coverage_fraction column
df$coverage_fraction <- NULL

You can learn more about this package here: https://github.com/isciences/exactextractr