0
votes

I'm trying to find the mean daily temperature for counties in South Dakota from raster grids ('bil' files) found at http://prism.oregonstate.edu/. I am getting county boundaries from the 'maps' package.

library(maps)
library(raster)
sd_counties <- map('county','south dakota')
sd_raster <- raster('file_path')

How do I extract the grid cells within each county? I think I need to turn each county into it's own polygon to do this, but how? Then, I should be able to do something like the following. Any help would be greatly appreciated.

values <- extract(raster, list of polygons)  
polygon_means <- unlist(lapply(values, FUN=mean))
1

1 Answers

1
votes

I'm not familiar with the maps package or the map function, but it looks like it's solely for visualization, rather than geospatial operations.

While there might be a way to convert the map object to actual polygons, here's an easy way sing raster's getData function that works:

library(raster)

usa_adm2 <- getData(country='USA',level=2)

sd_counties <- usa_adm2[grepl('South Dakota',usa_adm2$NAME_1),]

plot(sd_counties)

enter image description here

Now you can extract pixels for each county using extract(r,sd_counties), where r is your desired raster.

Note, that depending on the number of pixels (and layers) you need to extract, that can take some time.