[R-sig-Geo] over() in sp package
Kenny Bell
kmb56 at berkeley.edu
Tue Apr 14 18:58:01 CEST 2015
Thanks so much for all this. Yes, I'm after the area weighted averages. The
eventual application is to aggregate a 2.5 arcminute temperature grid up to
a zip-5 polygon, many of which overlap only 2-3 grid cells.
On Mon, Apr 13, 2015 at 10:47 PM, Edzer Pebesma <
edzer.pebesma at uni-muenster.de> wrote:
>
>
> On 04/14/2015 01:27 AM, Kenny Bell wrote:
> > Hi,
> >
> > When I call:
> >
> > over(spPolyDf, spGridDf)
> >
> > The function seems to be returning the spGridDf$band1 value for the first
> > grid cell that is contained in each polygon.
>
> This is documented; you can use returnList = TRUE to get the list of
> values for each polygon. Since you want to aggregate these anyway, you
> might as well use
>
> aggregate(spGridDF, spPolyDf, mean)
>
> this is gives a simple average of pixels whose centers fall inside the
> polygons. Do you need a weighted average, weighted by how much of the
> pixel overlaps with each polygon?
>
> >
> > What I would like is an overlap-area weighted average (or some provided
> > function like sum) of all the grid cells contained (or partially
> contained)
> > in each polygon element of spPolyDf.
> >
> > I have code that does this by eventually using joinPolys() from
> PBSMapping,
> > but it's slow and would be much cleaner to have a sp-only solution. Is
> > there any way to do this? Some example code is below:
> >
> > library(sp)
> > library(rgdal)
> >
> > # Download a USA shapefile
> > tmpFile <- tempfile()
> > download.file("
> > http://www2.census.gov/geo/tiger/GENZ2013/cb_2013_us_state_20m.zip",
> > destfile = tmpFile)
> > unzipped <- unzip(zipfile = tmpFile, exdir = getwd())
> >
> > # This is the actual reading stage
> > states <- readOGR(getwd(), layer="cb_2013_us_state_20m")
> >
> > # Remove all the files on disk
> > file.remove(tmpFile)
> > file.remove(unzipped)
> >
> > tmpFile <- tempfile()
> >
> > url <- paste0("http://www.prism.oregonstate.edu/fetchData.php?type",
> > "=bil&kind=recent&elem=tmax&range=daily&temporal=",
> > "20140416")
> >
> > # Download the file to fileName
> > download.file(url = url, mode = "wb", destfile = tmpFile)
> >
> > # unzip the file to the working directory.
> > unzippedFiles <- unzip(zipfile = tmpFile, overwrite = TRUE)
> >
> > # read into R
> > tMaxGrid <- readGDAL(fname = unzippedFiles[1])
> >
> > # delete the file
> > file.remove(tmpFile)
> > file.remove(unzippedFiles)
> >
> > # This doesn't work correctly
> > states$tmax20140416 <- over(states, tMaxGrid)
>
> this would assign a data.frame to what you want to be a vector (in a
> data.frame), so it would have to be:
>
> states$tmax20140416 <- over(states, tMaxGrid)[[1]]
>
> I'll try to modify sp so that it'd warn you in this case.
>
>
> --
> Edzer Pebesma
> Institute for Geoinformatics (ifgi), University of Münster,
> Heisenbergstraße 2, 48149 Münster, Germany; +49 251 83 33081
> Journal of Statistical Software: http://www.jstatsoft.org/
> Computers & Geosciences: http://elsevier.com/locate/cageo/
> Spatial Statistics Society http://www.spatialstatistics.info
>
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>
--
Kendon Bell
Email: kmb56 at berkeley.edu
Phone: (510) 612-3375
Ph.D. Student
Department of Agricultural & Resource Economics
University of California, Berkeley
Graduate Student Researcher
Energy Biosciences Institute
http://www.energybiosciencesinstitute.org/
[[alternative HTML version deleted]]
More information about the R-sig-Geo
mailing list