22

I have a set of raster images which represent a specific month over the years, and I want to make a timeline of the average values of an area using a shapefile.

How do I extract the values from the rasters and import them in R in a manner that I can use them?

3 Answers 3

28

Here is example code. It is fairly straight forward to adapt this code to work in a loop for processing all of your rasters. If your rasters share a common extent and resolution you can create a raster stack and loop through the bands in the stack. To create a vector containing all rasters in a directory, in a specific format, you can use "list.files" and then pass this vector to stack.

Example:

rlist=list.files(getwd(), pattern="img$", full.names=TRUE) 

r <- stack(rlist)   


    # Add required libraries
    require(raster)
    require(sp)
    require(rgdal)

    # Set working directory, raster, in and out shapefiles
    setwd("C:/test")
    inshp="MyPolys"
    outshp="PolyMeans"
    rdata <- "Year2012.img"

    # Read polygon feature class shapefile
    sdata <- readOGR(dsn=getwd(), layer=inshp)

    # Read raster
    r <- raster(rdata)

    # Extract raster values to list object
    r.vals <- extract(r, sdata)

    # Use list apply to calculate mean for each polygon
    r.mean <- lapply(r.vals, FUN=mean)

    # Join mean values to polygon data
    sdata@data <- data.frame(sdata@data, m2012=r.mean)

    # Write results
    writeOGR(sdata, getwd(), outshp, driver="ESRI Shapefile", check_exists=TRUE, 
         overwrite_layer=TRUE)
3
  • 1
    +1 - excellent solution, well laid out and complete answer.
    – Simbamangu
    Jul 10, 2012 at 4:10
  • 1
    If I have a multipolygonfile: How does the data.frame(sdata@data, m2012=r.mean) know whch polygon to assign which value to?
    – four-eyes
    Mar 10, 2016 at 13:11
  • 2
    to get r.mean added correctly to sdata I had to first unlist r.mean: r.mean <- unlist(lapply(r.vals, FUN=mean))
    – cmbarbu
    Nov 28, 2018 at 13:52
7

Read the shapefile into a SpatialPolygonsDataFrame (readOGR function from package rgdal)

Read the raster into a Raster object (raster function from package raster)

Use extract(raster, spdf) to get the grid cells under each polygon. Then run mean on them.

Repeat over your set of raster images...

1
  • how to write the shapefile out (with the values extracted from the raster image)?
    – four-eyes
    Mar 10, 2016 at 13:11
3

Both answer above are great, I wanted to point out that the raster::extract() function does allow you to call functions within it.

r1 <- raster
p1 <- polygon feature

p1$mean_val <- raster::extract(x = r1, y = p1, FUN=mean) 

Makes for a nice one liner for adding your summary values back to your polygon feature.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.