extract: Extract values from Raster objects

Description

Extract data from a Raster object using cell numbers or the locations of other spatial data (i.e., a spatial query). You can use coordinates (points), lines, polygons or an Extent (rectangle) object.

Usage

extract(x, y, ...)

Arguments

Raster* object

A vector (representing cell numbers); or a SpatialPoints object or a two-column data.frame or matrix (representing points); or a SpatialPolygons, S

...

Additional arguments, see under Details

Value

A vector for RasterLayer objects, and a matrix for a RasterStack or RasterBrick object. A list if y is a SpatialPolygons* or SpatialLines* object or if a buffer argument is used (but not a fun argument). focal values are returned as a n column matrix. The first column has the column numbers, the remaining has the focal values for each layer.

Details

Below are additional arguments that can be used, depending on the type of objects used. = all objects = 1) fun. Function to summarize the values (e.g. mean). The function should take a single numeric vector as argument and return a single values (e.g. mean, min or max), and accept a na.rm argument. Thus, standard R functions not including an na.rm argument must be wrapped as in this example: fun=function(x,...)length(x) If y represents points supplying a fun argument is only useful when a buffer is used (see below). 2) na.rm. Only useful when fun is supplied. If na.rm=TRUE (the default value), NA values are removed before fun is applied. This argument may be ignored if the function used has a ... argument and ignres an additional na.rm argument. = points = If y represents points, extract returns the values of a Raster* object for the cells in which a set of points fall. Additional arguments that can be suplied to this function if y represents points: 1) method. If method='simple' (the default), values for the cell a point falls in are returned. The alternative is method='bilinear', in which case the returned values are interpolated from the values of the four nearest raster cells. 2) buffer. The radius of a buffer around each point for which to extract cell values. If the distance between the sampling point and the center of a cell is less than or equal to the buffer, the cell is included. The buffer can be specified as a single value, or as a vector of the length of the number of points. If the data are not projected (latitude/longitude), the unit should be meters. Otherwise it should be in map-units (typically also meters). 3) cellnumbers. Logical. If cellnumbers=TRUE and buffer is not NULL, cellnumbers will also be returned. cellnumbers is set to FALSE when an argument fun is used. = lines = df. Logical. If df=TRUE, results will be returned as a data.frame, rather than as a list. And, if no fun argument is supplied: cellnumbers. Logical. If cellnumbers=TRUE is used, cell-numbers will also be returned. = polygons = If y represents polygons, the extract method returns the values of the cells of a Raster* object that is covered by a polygon. A cell is covered if its center is inside the polygon (but see the weights option for considering partly covered cells; and argument small for getting values for small polygons anyway). if y represents polygons are: 1) small. Logical. Default is FALSE. If TRUE a value is also returned for relatively small polygons (e.g. those smaller than a single cell of the Raster* object), or polygons with an odd shape, for which otherwise no values are returned because they do not cover any Raster* object's cell centers. In some cases, you could use alternatively use the centroids of such polygons, for example using extract(x, coordinates(y)) or extract(x, coordinates(y), method='bilinear'). 2) df. Logical. If df=TRUE, results will be returned as a data.frame, rather than as a list. Additional arguments that can be suplied to this function, if no fun argument is supplied: 3) weights. If TRUE, the function returns, for each polygon, a matrix with the cell values and the approximate fraction of each cell that is covered by the polygon(rounded to 1/100). The weights can be used for averaging; see examples. This option can be useful if the polygons are small relative to the cells size of the Raster* object. 4) cellnumbers. Logical. If cellnumbers=TRUE, cell-numbers will also be returned. = focal values = If y is missing, you can supply two arguments to get 'focal' values (i.e. the values in a neighborhood around each cell) 1) row. row number (this is a required argument). 2) ngb. The neighborhood to consider. Either a single integer or a vector of two integers. Default is 3. See focal = all multi-layer objects = If x is a RasterStack or RasterBrick object, extract accepts, in addition to the arguments listed above, these two additional arguments: 1)layer. Integer. First layer for which you want values 2) nl. Integer. Number of layers for which you want values

Examples

Run this code

r <- raster(ncol=36, nrow=18)
r[] <- 1:ncell(r)

###############################
# extract values by cell number
###############################
extract(r, c(1:2, 10, 100))
s <- stack(r, sqrt(r), r/r)
extract(s, c(1, 10, 100), layer=2, n=2)

###############################
# extract values with points
###############################
xy <- cbind(-50, seq(-80, 80, by=20))
extract(r, xy)

sp <- SpatialPoints(xy)
extract(r, sp, method='bilinear')

# examples with a buffer
extract(r, xy[1:3,], buffer=1000000)
extract(r, xy[1:3,], buffer=1000000, fun=mean)

## illustrating the varying size of a buffer (expressed in meters) on a longitude/latitude raster
 z <- extract(r, xy, buffer=1000000)
 s <- raster(r)
 for (i in 1:length(z)) { s[z[[i]]] <- i }
## compare with raster that is not longitude/latitude
 projection(r) <- "+proj=utm +zone=17" 
 xy[,1] <- 50
 z <- extract(r, xy, buffer=8)
 for (i in 1:length(z)) { s[z[[i]]] <- i }
 plot(s)
# library(maptools)
# data(wrld_simpl)
# plot(wrld_simpl, add=TRUE)

###############################
# extract values with lines
###############################

cds1 <- rbind(c(-50,0), c(0,60), c(40,5), c(15,-45), c(-10,-25))
cds2 <- rbind(c(80,20), c(140,60), c(160,0), c(140,-55))
lines <- SpatialLines(list(Lines(list(Line(cds1)), "1"), Lines(list(Line(cds2)), "2") ))

extract(r, lines)

###############################
# extract values with polygons
###############################
cds1 <- rbind(c(-180,-20), c(-160,5), c(-60, 0), c(-160,-60), c(-180,-20))
cds2 <- rbind(c(80,0), c(100,60), c(120,0), c(120,-55), c(80,0))
polys <- SpatialPolygons(list(Polygons(list(Polygon(cds1)), 1), Polygons(list(Polygon(cds2)), 2)))

#plot(r)
#plot(polys, add=TRUE)
v <- extract(r, polys)
v
# mean for each polygon
unlist(lapply(v, function(x) if (!is.null(x)) mean(x, na.rm=TRUE) else NA ))

# v <- extract(r, polys, cellnumbers=TRUE)

# weighted mean
# v <- extract(r, polys, weights=TRUE, fun=mean)
# equivalent to:
# v <- extract(r, polys, weights=TRUE)
# sapply(v, function(x) if (!is.null(x)) {sum(apply(x, 1, prod)) / sum(x[,2])} else NA  )



###############################
# extract values with an extent
###############################
e <- extent(150,170,-60,-40)
extract(r, e)
#plot(r)
#plot(e, add=T)


###############################
# focal values
###############################
f <- extract(r, row=5, ngb=3)
f[1:15, ]