Extract raster returns NA, when a polygon is inside a single pixel

[hriihimaki]

I noticed this issue when I was extracting raster data for 6 m radius plots with raster resolution from 1 to 30 meters. At the coarser resolution some of the polygons do not coincide with any pixel centers. Thus, the function (e.g. mean or median) returns NA:

# Reproducible example:   
library(velox)
library(raster)

## Make VeloxRaster
mat <- matrix(1:100, 10, 10)
extent <- c(0,1,0,1)
vx <- velox(mat, extent=extent, res=c(0.1,0.1), crs="+proj=longlat +datum=WGS84 +no_defs")

## Make SpatialPolygonsDataFrame
library(sp)
library(rgeos)
set.seed(0)
coords <- cbind(runif(10, extent[1], extent[2]), runif(10, extent[3], extent[4]))
sp <- SpatialPoints(coords)

# Default example
# from https://cran.r-project.org/web/packages/velox/README.html
spol_norm <- gBuffer(sp, width=0.2, byid=TRUE)
spdf_norm <- SpatialPolygonsDataFrame(spol_norm, data.frame(id=1:length(spol_norm)), FALSE)

# Smaller buffer
spol_small<- gBuffer(sp, width=0.05, byid=TRUE)
spdf_small <- SpatialPolygonsDataFrame(spol_small, data.frame(id=1:length(spol_small)), FALSE)

plot(spdf_norm); par(new=F)
plot(spdf_small)

## Extract values and calculate mean, see results
(ex.mat.norm <- vx$extract(spdf_norm, fun=median))
(ex.mat.small <- vx$extract(spdf_small, fun=median)) # -> 3 NA's


# Check problematic polygons
r <- raster(mat)
plot(r, main="Problematic polygons")
lines(spdf_small) 
lines(spdf_small[5:7,], col="red") # Plots NAs in red

[hunzikp]

This problem can now be addressed by setting argument small = TRUE in VeloxRaster_extract; in this case, small (or oddly shaped) polygons are intersected with the full grid cells, instead of only the grid cell centers.

Applied to your example:

# Reproducible example:   
library(velox)
library(raster)

## Make VeloxRaster
mat <- matrix(1:100, 10, 10)
extent <- c(0,1,0,1)
vx <- velox(mat, extent=extent, res=c(0.1,0.1), crs="+proj=longlat +datum=WGS84 +no_defs")

## Make SpatialPolygonsDataFrame
library(sp)
library(rgeos)
set.seed(0)
coords <- cbind(runif(10, extent[1], extent[2]), runif(10, extent[3], extent[4]))
sp <- SpatialPoints(coords)

# Default example
# from https://cran.r-project.org/web/packages/velox/README.html
spol_norm <- gBuffer(sp, width=0.2, byid=TRUE)
spdf_norm <- SpatialPolygonsDataFrame(spol_norm, data.frame(id=1:length(spol_norm)), FALSE)

# Smaller buffer
spol_small<- gBuffer(sp, width=0.05, byid=TRUE)
spdf_small <- SpatialPolygonsDataFrame(spol_small, data.frame(id=1:length(spol_small)), FALSE)

plot(spdf_norm); par(new=F)
plot(spdf_small)

## Extract values and calculate mean, see results
(ex.mat.norm <- vx$extract(spdf_norm, fun=median))
(ex.mat.small <- vx$extract(spdf_small, fun=median, small = TRUE)) # -> No NA values anymore.