# Instantiating spatial polygon without using a shapefile in R

So, the usual way we read a shapefile in R is via the maptools package, like this:

However, I have a use case whereby I don't have a shapefile.shp but instead I have a series of polygon coordinates

16.484375 59.736328125,17.4951171875 55.1220703125,24.74609375 55.0341796875,22.5927734375 61.142578125,16.484375 59.736328125

and its corresponding projection

coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0

How do I "instantiate" sfdata (which will be a "polygon object") directly from this data? (without going in a roundabout way of creating a shapefile with these data and then reading from the newly created shapefile)

First get the coordinates into a 2-column matrix:

> xym
[,1]     [,2]
[1,] 16.48438 59.73633
[2,] 17.49512 55.12207
[3,] 24.74609 55.03418
[4,] 22.59277 61.14258
[5,] 16.48438 59.73633

Then create a Polygon, wrap that into a Polygons object, then wrap that into a SpatialPolygons object:

> library(sp)
> p = Polygon(xym)
> ps = Polygons(list(p),1)
> sps = SpatialPolygons(list(ps))

The reason for this level of complexity is that a Polygon is a simple ring, a Polygons object can be several rings with an ID (here set to 1) (so is like a single feature in a GIS) and a SpatialPolygons can have a CRS. Ooh, I should probably set it:

> proj4string(sps) = CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0")

If you want to turn it into a SpatialPolygonsDataFrame (which is what comes our of readShapeSpatial when the shapefile is polygons) then do:

> data = data.frame(f=99.9)
> spdf = SpatialPolygonsDataFrame(sps,data)
> spdf

giving this:

> summary(spdf)
Object of class SpatialPolygonsDataFrame
Coordinates:
min      max
x 16.48438 24.74609
y 55.03418 61.14258
Is projected: FALSE
proj4string :
[+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0]
Data attributes:
Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
99.9    99.9    99.9    99.9    99.9    99.9
• +1 Very nice, clear exposition. It's great to see the code broken up by explanations rather than being offered as a monolithic block! – whuber Dec 30 '11 at 14:56
• Excellent ... great to see how these objects are put together! Need to see more of the R help pages written clearly like this. – Simbamangu Dec 30 '11 at 17:26
• Its something I have to re-teach myself every time I want to do it, so I take any opportunity to teach other people! – Spacedman Dec 31 '11 at 13:13
• excellent... how would I go about adding multiple unique id (f) polygons to the data frame? – mga Jun 26 '13 at 15:39
• For this answer to have more general validity, could you show how to do it in case of multiple polygons? This is a bit tricky. – Tomas Oct 23 '14 at 10:36

To complete Spacedman's excellent answer for the case where your data would contain multiple polygons, here is some code using dplyr:

library(dplyr)
library(ggplot2)
library(sp)
## use data from ggplot2:::geom_polygon example:
positions <- data.frame(id = rep(factor(c("1.1", "2.1", "1.2", "2.2", "1.3", "2.3")), each = 4),
x = c(2, 1, 1.1, 2.2, 1, 0, 0.3, 1.1, 2.2, 1.1, 1.2, 2.5, 1.1, 0.3,
0.5, 1.2, 2.5, 1.2, 1.3, 2.7, 1.2, 0.5, 0.6, 1.3),
y = c(-0.5, 0, 1, 0.5, 0, 0.5, 1.5, 1, 0.5, 1, 2.1, 1.7, 1, 1.5,
2.2, 2.1, 1.7, 2.1, 3.2, 2.8, 2.1, 2.2, 3.3, 3.2)) %>% as.tbl

df_to_spp <- positions %>%
group_by(id) %>%
do(poly=select(., x, y) %>%Polygon()) %>%
rowwise() %>%
do(polys=Polygons(list(.\$poly),.\$id)) %>%
{SpatialPolygons(.\$polys)}

## plot it
plot(df_to_spp)

Just for fun, you can compare with the plot obtained with ggplot2 using the initial data-frame:

ggplot(positions) +
geom_polygon(aes(x=x, y=y, group=id), colour="black", fill=NA)

Note that the code above assumes that you have only one polyogn per id. If some ids had disjoint polygons, I guess one should add another column in the dataset, first group_by the sub-id, then use group_by(upper-id) instead of rowwise

Same code using the purrr::map function:

df_to_spp <- positions %>%
nest(-id) %>%
mutate(Poly=purrr::map(data, ~select(., x, y)  %>% Polygon()),
polys=map2(Poly, id, ~Polygons(list(.x),.y))) %>%
{SpatialPolygons(.\$polys)}