2

I have a dataframe with around 40 points. They are grouped by an ID variable (their name), which four points in each group. I'm trying to turn these points into polygons, but I've had a lot of trouble with different methods in both sf and sp — currently I have them as a large sf dataframe, but I can't seem to find any function that would use a grouping variable to create polygons.

I've tried using the st_cast method, the general idea being

polypoints <- polypoints %>% group_by("ID") %>% st_cast("MULTIPOINT") %>% st_cast("MULTILINESTRING") %>% st_cast("MULTIPOLYGON")

but the step from multipoint to multiline string leaves the geography empty — what are the steps for a process like this?

  • is MULITPOLYGON a typo in the code or just here? – Spacedman Aug 19 at 19:51
  • please try and give us a reproducible example with sample data? – Spacedman Aug 19 at 19:52
  • yes the "MULITPOLYGON" is just a type, I just fixed it — I'll update with some sample data shortly – cschwab98 Aug 19 at 20:02
  • Never mind think I've done it all for you see below – Spacedman Aug 19 at 20:03
2

Make a reproducible example:

> set.seed(999)
> xy = data.frame(x=runif(24), y=runif(24))
> xy$ID = rep(c(1:6), rep(4,6))
> head(xy)
           x         y ID
1 0.38907138 0.8260703  1
2 0.58306072 0.8195141  1
3 0.09466569 0.5684927  1
4 0.85263123 0.6196068  1
5 0.78674676 0.8308805  2
6 0.11934226 0.4588336  2

Make the data frame into an sf data frame:

> xys = st_as_sf(xy, coords=c("x","y"))

Then aggregate by ID, combine the points and cast to POLYGON, turn the whole thing into an sf data frame. In a one-liner:

> polys = st_sf(
  aggregate(
    xys$geometry,
    list(xys$ID),
    function(g){
       st_cast(st_combine(g),"POLYGON")
    }
   ))

> polys
Simple feature collection with 6 features and 1 field
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0.03014558 ymin: 0.01902308 xmax: 0.9875201 ymax: 0.8690149
epsg (SRID):    NA
proj4string:    NA
  Group.1                       geometry
1       1 POLYGON ((0.3890714 0.82607...
2       2 POLYGON ((0.7867468 0.83088...
3       3 POLYGON ((0.3907724 0.52724...
4       4 POLYGON ((0.03014558 0.8162...
5       5 POLYGON ((0.1665847 0.10961...
6       6 POLYGON ((0.9074913 0.35951...

Looks pretty awful but that's random data for you...

enter image description here

  • thank you so much, this is extremely helpful (both in solving the problem and in teaching me new ways to think about gis in R) - I have a few quick questions though. Where is the xys dataframe coming from? Also, trying to use this on my own data returns: Error in st_sfc(NextMethod(), crs = st_crs(x), precision = st_precision(x)): is.numeric(crs) || is.character(crs) || inherits(crs, "crs") is not TRUE . Not sure what is causing this, I tried including st_crs in the polys line but that didn't help either. – cschwab98 Aug 19 at 20:22
  • Just to add on to the error, using my own data won't even display the object using view("polys"), and says r error 4 (R code execution error) despite compiling properly – cschwab98 Aug 19 at 20:26
  • Oops missing line - xys is an sf data frame created from xy - edited in. – Spacedman Aug 19 at 20:29
  • yes, figured that one out! still stuck on crs issue, which doesn't happen in sample issue. Not sure what's happening, tried not including crs at st_as_sf step but that doesn't help – cschwab98 Aug 19 at 20:35
  • Can't duplicate. If I st_crs(xys)=4326 and then run it it works, I can even add ,crs=4326) in the wrapping st_sf( call.... – Spacedman Aug 19 at 20:43

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.