# Grouping polygons based on clusters

I have a polygon `sf` and all I want is to get the polygons grouped based on the clusters they belong to - ideally, the cluster grouping should be added as one more attribute to each polygon. The definition of cluster here is any polygon that shares an edge or a vertex with any neighbouring polygon will be in the same cluster.

My attemp with using `poly2ng` is as follows

``````tiles_cluster <- poly2nb(tiles)
coords <- st_coordinates(st_centroid(st_geometry(tiles)))
plot(st_geometry(tiles), border="grey")
plot(tiles_cluster, coords, col="grey", add = TRUE)
``````

But I don't know how to extract the clusters from the output.

I would also like to know if I can cluster them using a threshold distance meaning even if the polygons are discountinous if they are within the threshold distance they should be clustered together.

``````tiles

structure(list(x = c(15L, 14L, 15L, 16L, 14L, 15L, 16L, 17L,
13L, 14L, 15L, 16L, 17L, 18L, 12L, 13L, 14L, 15L, 16L, 17L, 11L,
12L, 13L, 14L, 15L, 16L, 11L, 13L, 14L, 15L, 16L, 13L, 14L, 15L,
13L, 14L, 11L, 8L, 9L, 10L, 11L, 12L, 7L, 8L, 9L, 10L, 11L, 12L,
13L, 6L, 7L, 11L, 12L, 13L, 7L, 8L, 12L, 13L, 8L, 14L, 8L, 15L,
8L, 16L, 7L, 15L, 14L, 15L, 16L, 14L, 15L, 16L, 17L, 13L, 14L
), y = c(2L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L,
5L, 6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 7L, 7L, 8L, 8L, 8L,
8L, 8L, 9L, 9L, 9L, 10L, 10L, 11L, 12L, 12L, 12L, 12L, 12L, 13L,
13L, 13L, 13L, 13L, 13L, 13L, 14L, 14L, 14L, 14L, 14L, 15L, 15L,
15L, 15L, 16L, 16L, 17L, 17L, 18L, 18L, 20L, 2L, 3L, 3L, 3L,
4L, 4L, 4L, 4L, 5L, 5L), haz = c(0.226500004529953, 0.226899996399879,
0.268900007009506, 0.342599987983704, 0.235200002789497, 0.342099994421005,
0.990299999713898, 0.580600023269653, 0.228200003504753, 0.324699997901917,
1.32990002632141, 1.70360004901886, 0.628899991512299, 0.199000000953674,
0.153699994087219, 0.25450000166893, 0.156499996781349, 0.79720002412796,
0.8367999792099, 0.196700006723404, 0.109300002455711, 0.209600001573563,
0.174199998378754, 0.20440000295639, 0.359400004148483, 0.286500006914139,
0.143500000238419, 0.233400002121925, 0.32480001449585, 0.359800010919571,
0.164100006222725, 0.233500003814697, 0.347799986600876, 0.3125,
0.236000001430511, 0.267500013113022, 0.107600003480911, 0.267500013113022,
0.602299988269806, 0.527700006961823, 0.329199999570847, 0.200399994850159,
0.205300003290176, 0.110699996352196, 0.5242999792099, 1.16700005531311,
0.865100026130676, 0.370400011539459, 0.184400007128716, 0.211600005626678,
0.214300006628036, 0.620500028133392, 0.561200022697449, 0.364499986171722,
0.145699992775917, 0.177100002765656, 0.207599997520447, 0.213899999856949,
0.196600005030632, 0.439900010824203, 0.288800001144409, 0.311500012874603,
0.206900000572205, 0.10419999808073, 0.137799993157387, 0.226500004529953,
0.226899996399879, 0.268900007009506, 0.342599987983704, 0.235200002789497,
0.342099994421005, 0.990299999713898, 0.580600023269653, 0.228200003504753,
0.324699997901917), geom = structure(list(structure(list(structure(c(384875,
384880, 384880, 384875, 384875, 205270, 205270, 205265, 205265,
205270), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384870, 384875, 384875, 384870,
384870, 205265, 205265, 205260, 205260, 205265), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384875, 384880, 384880, 384875, 384875, 205265,
205265, 205260, 205260, 205265), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384880, 384885,
384885, 384880, 384880, 205265, 205265, 205260, 205260, 205265
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384870, 384875, 384875, 384870,
384870, 205260, 205260, 205255, 205255, 205260), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384875, 384880, 384880, 384875, 384875, 205260,
205260, 205255, 205255, 205260), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384880, 384885,
384885, 384880, 384880, 205260, 205260, 205255, 205255, 205260
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384885, 384890, 384890, 384885,
384885, 205260, 205260, 205255, 205255, 205260), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384865, 384870, 384870, 384865, 384865, 205255,
205255, 205250, 205250, 205255), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384870, 384875,
384875, 384870, 384870, 205255, 205255, 205250, 205250, 205255
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384875, 384880, 384880, 384875,
384875, 205255, 205255, 205250, 205250, 205255), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384880, 384885, 384885, 384880, 384880, 205255,
205255, 205250, 205250, 205255), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384885, 384890,
384890, 384885, 384885, 205255, 205255, 205250, 205250, 205255
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384890, 384895, 384895, 384890,
384890, 205255, 205255, 205250, 205250, 205255), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384860, 384865, 384865, 384860, 384860, 205250,
205250, 205245, 205245, 205250), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384865, 384870,
384870, 384865, 384865, 205250, 205250, 205245, 205245, 205250
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384870, 384875, 384875, 384870,
384870, 205250, 205250, 205245, 205245, 205250), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384875, 384880, 384880, 384875, 384875, 205250,
205250, 205245, 205245, 205250), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384880, 384885,
384885, 384880, 384880, 205250, 205250, 205245, 205245, 205250
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384885, 384890, 384890, 384885,
384885, 205250, 205250, 205245, 205245, 205250), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384855, 384860, 384860, 384855, 384855, 205245,
205245, 205240, 205240, 205245), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384860, 384865,
384865, 384860, 384860, 205245, 205245, 205240, 205240, 205245
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384865, 384870, 384870, 384865,
384865, 205245, 205245, 205240, 205240, 205245), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384870, 384875, 384875, 384870, 384870, 205245,
205245, 205240, 205240, 205245), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384875, 384880,
384880, 384875, 384875, 205245, 205245, 205240, 205240, 205245
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384880, 384885, 384885, 384880,
384880, 205245, 205245, 205240, 205240, 205245), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384855, 384860, 384860, 384855, 384855, 205240,
205240, 205235, 205235, 205240), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384865, 384870,
384870, 384865, 384865, 205240, 205240, 205235, 205235, 205240
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384870, 384875, 384875, 384870,
384870, 205240, 205240, 205235, 205235, 205240), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384875, 384880, 384880, 384875, 384875, 205240,
205240, 205235, 205235, 205240), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384880, 384885,
384885, 384880, 384880, 205240, 205240, 205235, 205235, 205240
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384865, 384870, 384870, 384865,
384865, 205235, 205235, 205230, 205230, 205235), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384870, 384875, 384875, 384870, 384870, 205235,
205235, 205230, 205230, 205235), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384875, 384880,
384880, 384875, 384875, 205235, 205235, 205230, 205230, 205235
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384865, 384870, 384870, 384865,
384865, 205230, 205230, 205225, 205225, 205230), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384870, 384875, 384875, 384870, 384870, 205230,
205230, 205225, 205225, 205230), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384855, 384860,
384860, 384855, 384855, 205225, 205225, 205220, 205220, 205225
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384840, 384845, 384845, 384840,
384840, 205220, 205220, 205215, 205215, 205220), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384845, 384850, 384850, 384845, 384845, 205220,
205220, 205215, 205215, 205220), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384850, 384855,
384855, 384850, 384850, 205220, 205220, 205215, 205215, 205220
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384855, 384860, 384860, 384855,
384855, 205220, 205220, 205215, 205215, 205220), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384860, 384865, 384865, 384860, 384860, 205220,
205220, 205215, 205215, 205220), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384835, 384840,
384840, 384835, 384835, 205215, 205215, 205210, 205210, 205215
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384840, 384845, 384845, 384840,
384840, 205215, 205215, 205210, 205210, 205215), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384845, 384850, 384850, 384845, 384845, 205215,
205215, 205210, 205210, 205215), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384850, 384855,
384855, 384850, 384850, 205215, 205215, 205210, 205210, 205215
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384855, 384860, 384860, 384855,
384855, 205215, 205215, 205210, 205210, 205215), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384860, 384865, 384865, 384860, 384860, 205215,
205215, 205210, 205210, 205215), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384865, 384870,
384870, 384865, 384865, 205215, 205215, 205210, 205210, 205215
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384830, 384835, 384835, 384830,
384830, 205210, 205210, 205205, 205205, 205210), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384835, 384840, 384840, 384835, 384835, 205210,
205210, 205205, 205205, 205210), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384855, 384860,
384860, 384855, 384855, 205210, 205210, 205205, 205205, 205210
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384860, 384865, 384865, 384860,
384860, 205210, 205210, 205205, 205205, 205210), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384865, 384870, 384870, 384865, 384865, 205210,
205210, 205205, 205205, 205210), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384835, 384840,
384840, 384835, 384835, 205205, 205205, 205200, 205200, 205205
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384840, 384845, 384845, 384840,
384840, 205205, 205205, 205200, 205200, 205205), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384860, 384865, 384865, 384860, 384860, 205205,
205205, 205200, 205200, 205205), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384865, 384870,
384870, 384865, 384865, 205205, 205205, 205200, 205200, 205205
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384840, 384845, 384845, 384840,
384840, 205200, 205200, 205195, 205195, 205200), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384870, 384875, 384875, 384870, 384870, 205200,
205200, 205195, 205195, 205200), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384840, 384845,
384845, 384840, 384840, 205195, 205195, 205190, 205190, 205195
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384875, 384880, 384880, 384875,
384875, 205195, 205195, 205190, 205190, 205195), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384840, 384845, 384845, 384840, 384840, 205190,
205190, 205185, 205185, 205190), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384880, 384885,
384885, 384880, 384880, 205190, 205190, 205185, 205185, 205190
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384835, 384840, 384840, 384835,
384835, 205180, 205180, 205175, 205175, 205180), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384875, 384880, 384880, 384875, 384875, 205270,
205270, 205265, 205265, 205270), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384870, 384875,
384875, 384870, 384870, 205265, 205265, 205260, 205260, 205265
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384875, 384880, 384880, 384875,
384875, 205265, 205265, 205260, 205260, 205265), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384880, 384885, 384885, 384880, 384880, 205265,
205265, 205260, 205260, 205265), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384870, 384875,
384875, 384870, 384870, 205260, 205260, 205255, 205255, 205260
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384875, 384880, 384880, 384875,
384875, 205260, 205260, 205255, 205255, 205260), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384880, 384885, 384885, 384880, 384880, 205260,
205260, 205255, 205255, 205260), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg")), structure(list(structure(c(384885, 384890,
384890, 384885, 384885, 205260, 205260, 205255, 205255, 205260
), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")),
structure(list(structure(c(384865, 384870, 384870, 384865,
384865, 205255, 205255, 205250, 205250, 205255), dim = c(5L,
2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
structure(c(384870, 384875, 384875, 384870, 384870, 205255,
205255, 205250, 205250, 205255), dim = c(5L, 2L))), class = c("XY",
"POLYGON", "sfg"))), class = c("sfc_POLYGON", "sfc"), precision = 0, bbox = structure(c(xmin = 384830,
ymin = 205175, xmax = 384895, ymax = 205270), class = "bbox"), crs = structure(list(
input = "EPSG:27700", wkt = "PROJCRS[\"OSGB36 / British National Grid\",\n    BASEGEOGCRS[\"OSGB36\",\n        DATUM[\"Ordnance Survey of Great Britain 1936\",\n            ELLIPSOID[\"Airy 1830\",6377563.396,299.3249646,\n                LENGTHUNIT[\"metre\",1]]],\n        PRIMEM[\"Greenwich\",0,\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n        ID[\"EPSG\",4277]],\n    CONVERSION[\"British National Grid\",\n        METHOD[\"Transverse Mercator\",\n            ID[\"EPSG\",9807]],\n        PARAMETER[\"Latitude of natural origin\",49,\n            ANGLEUNIT[\"degree\",0.0174532925199433],\n            ID[\"EPSG\",8801]],\n        PARAMETER[\"Longitude of natural origin\",-2,\n            ANGLEUNIT[\"degree\",0.0174532925199433],\n            ID[\"EPSG\",8802]],\n        PARAMETER[\"Scale factor at natural origin\",0.9996012717,\n            SCALEUNIT[\"unity\",1],\n            ID[\"EPSG\",8805]],\n        PARAMETER[\"False easting\",400000,\n            LENGTHUNIT[\"metre\",1],\n            ID[\"EPSG\",8806]],\n        PARAMETER[\"False northing\",-100000,\n            LENGTHUNIT[\"metre\",1],\n            ID[\"EPSG\",8807]]],\n    CS[Cartesian,2],\n        AXIS[\"(E)\",east,\n            ORDER[1],\n            LENGTHUNIT[\"metre\",1]],\n        AXIS[\"(N)\",north,\n            ORDER[2],\n            LENGTHUNIT[\"metre\",1]],\n    USAGE[\n        SCOPE[\"Engineering survey, topographic mapping.\"],\n        AREA[\"United Kingdom (UK) - offshore to boundary of UKCS within 49°45'N to 61°N and 9°W to 2°E; onshore Great Britain (England, Wales and Scotland). Isle of Man onshore.\"],\n        BBOX[49.75,-9,61.01,2.01]],\n    ID[\"EPSG\",27700]]"), class = "crs"), n_empty = 0L)), row.names = c(NA,
-75L), class = c("sf", "data.frame"), sf_column = "geom", agr = structure(c(x = NA_integer_,
y = NA_integer_, haz = NA_integer_), class = "factor", levels = c("constant",
"aggregate", "identity")))
``````
• Dissolve will convert them as single polygon per cluster I think. What I want is more complicated than that. Commented Dec 14, 2022 at 13:10

Try this:

``````require(sf)

geometries <- st_cast(st_union(st_buffer(shape, 0.1)), "POLYGON") #Buffer to make the corners touch, union to dissolve adjacent borders
dissolved <- st_sf(geometries) #Create a sf object from the geometries
dissolved\$clusterID=1:length(geometries) #Add cluster id to each row
result <- st_join(shape, dissolved) #Join cluster id to the original polygons

st_write(result, "/home/bera/Desktop/GIStest/grid_R_dissolved.shp")
``````

I'd like to show a different approach using {terra}, since this is easier to solve with raster data from my point of view:

``````library(terra)
#> terra 1.6.49

tiles <- vect(tiles)

# extract some attributes for next step
e <- ext(tiles)
ref <- crs(tiles)

# create the `y` attribute for `rasterize()`
r <- rast(extent = e, crs = ref, resolution = 5)

r2 <- rasterize(tiles, r, field = "haz")

# determine patches using Queen's case
patches <- patches(r2, directions = 8)

# inspect
plot(patches)
``````

Created on 2022-12-17 with reprex v2.0.2

• Thanks. This is a helpful answer, but I prefer BERA's answer as doesn't change the underlying data format Commented Feb 7, 2023 at 14:56