Dear Zivan,
Thank you for your email. I followed your scripts, only to find it a genius
solution.
Actually, only the first sentence is enough:
st_intersection(polygon) |> st_collection_extract('LINESTRING')
And it has brought me the interior boundary already.
However, I do not fully grasp the essence here. I know st_intersection(x, y) is
to find the shared part of x and y, but what is st_intersection(x)?
Conceptually it is "the shared part of x and x", but technically it brings me
more features than x itself (and I feel confused about the output). But
obviously you make the most of this feature.
Thank you once again!
Xiang
________________________________
From: R-sig-Geo <[email protected]> on behalf of Zivan Karaman
<[email protected]>
Sent: Thursday, August 7, 2025 03:07
To: Josiah Parry <[email protected]>
Cc: [email protected] <[email protected]>
Subject: Re: [R-sig-Geo] How to retrieve the interior boundary (points) of a
polygon sf object
Hi,
I'm not sure I understand what you want to achieve or why, but here's a try:
# Load required libraries
library(sf)
library(spData)
# Options to avoid messages from tigris
polygon <- st_geometry(us_states)
# Extract the shared borders only (where geometries touch more than once)
internal_lines <- st_intersection(polygon) |>
st_collection_extract("LINESTRING")
# Now remove the external boundary
# Step 1: Union all states into a single polygon (whole USA outline)
national_outline <- st_union(polygon) |> st_boundary()
# Step 2: Remove outer boundary from internal_lines
internal_only <- st_difference(internal_lines, national_outline)
# Plot result
plot(st_geometry(national_outline), col = "blue", lwd = 0.5, main =
"Interior Boundaries Only")
plot(st_geometry(internal_only), col = "red", lwd = 0.5, add = TRUE)
Best,
Zivan
On Wed, Aug 6, 2025 at 8:24 PM Josiah Parry <[email protected]> wrote:
>
> I'm not sure I fully understand the objective. But this is how I would do it
> with rsgeo—I'm not sure how to do it with sf
>
> library(sf)
> library(spData)
> library(rsgeo)
>
> union_geoms(as_rsgeo(us_states)) |>
> expand_geoms() |>
> .subset2(1) |>
> head(1) |>
> cast_geoms("multipoint") |>
> st_as_sfc() |>
> plot()
>
>
>
>
> On Wed, Aug 6, 2025 at 10:56 AM Andrea Gilardi - Unimib
> <[email protected]> wrote:
>>
>> Hi Xiang! I would use the following approach:
>>
>> library(tidyverse)
>> library(sf)
>> library(spData)
>> # You need the development version of sfnetworks:
>> # remotes::install_github("luukvdmeer/sfnetworks")
>> library(sfnetworks) # Exports st_duplicated
>>
>> polygon = st_geometry(us_states)
>> pts_matrix <- st_coordinates(polygon)
>> pts_sfc <- st_multipoint(pts_matrix[, c(1, 2)]) |> st_sfc(crs =
>> st_crs(polygon)) |> st_cast("POINT")
>> pts_dup <- st_duplicated(pts_sfc)
>> pts_interior <- pts_sfc[pts_dup]
>>
>> plot(polygon, reset = FALSE)
>> points(pts_interior, cex = 1, pch = 16, col = 2)
>>
>> You can also check the output here:
>> https://gist.github.com/agila5/2b5d9b7e1453c4d3f292e7aae5fdb79e
>>
>> Hope that's helpful.
>>
>> Andrea
>>
>> On 8/6/2025 7:27 PM, Dexter Locke wrote:
>> > library(tidyverse)
>> > library(sf)
>> > library(spData)
>> >
>> > polygon=st_geometry(us_states)
>> > plot(polygon)
>> >
>> > polygon |> st_as_sf() |> summarise() |> plot()
>> >
>> >
>> > On Wed, Aug 6, 2025 at 1:25 PM Xiang Ye <[email protected]> wrote:
>> >
>> >> Dear community,
>> >>
>> >> I am trapped by a seemingly easy task. I would like to write a short R
>> >> script to retrieve the interior boundary (points) of a polygon sf object,
>> >> like the state boundary of the Contiguous United States, but without the
>> >> outline.
>> >>
>> >> I tried several methods, with none fulfilling my objective. Below is the
>> >> code, with my observations in the comments:
>> >>
>> >> # Objective: Retrieve the internal boundary of a polygon data set
>> >>
>> >> library(tidyverse)
>> >> library(sf)
>> >> library(spData)
>> >>
>> >> polygon=st_geometry(us_states)
>> >> plot(polygon)
>> >>
>> >> # Method 1: Use st_difference()
>> >> st_union(polygon) -> border
>> >> plot(border)
>> >> st_difference(polygon, border) -> interior
>> >> plot(interior) # You can tell that the result is weird.
>> >>
>> >> # Method 2: Retrieve vertices as points, then use st_difference()
>> >> st_cast(polygon, 'POINT') -> points_polygon
>> >> points_polygon # 3585 features
>> >> st_cast(border, 'POINT') -> points_border
>> >> points_border # 1275 features
>> >> st_difference(points_polygon, points_border) -> points_interior # About 2
>> >> minutes to execute
>> >> points_interior # 4.57 million features, about 2 GB. why?
>> >>
>> >> # Methods 3: Use st_disjoint instead of st_difference()
>> >> polygon[border, op=st_disjoint] # Empty set
>> >> points_polygon[border, op=st_disjoint] -> points_interior2
>> >> points_interior2 # 132 features; plot() it and doesn't look right
>> >> points_polygon[points_border, op=st_disjoint] -> points_interior3
>> >> points_interior3 # 3585 features, meaning no points have been taken out
>> >>
>> >> # Methods 4: Retrieve vertices as coordinates, converting them into
>> >> tibbles, then anti-join the border coordinates
>> >> st_coordinates(polygon) %>% as_tibble -> coords_polygon
>> >> coords_polygon# A tibble: 3,585 × 5
>> >> st_coordinates(border) %>% as_tibble -> coords_border
>> >> coords_border # A tibble: 1,275 × 5
>> >> anti_join(coords_polygon, coords_border) -> coords_interior
>> >> coords_interior # A tibble: 3,585 × 5, meaning no points have been taken
>> >> out, and implying no coordinates in coords_border coincides with the
>> >> coordinates in coords_polygon. This is just not reasonable.
>> >>
>> >> Any comments on why I fail and/or how to do it correctly will be helpful.
>> >> Thank you in advance!
>> >>
>> >> 叶翔 YE, Xiang
>> >> THINKING SPATIALLY<http://www.linkedin.com/in/spatialyexiang>.
>> >> Ph.D. in Spatial Statistics
>> >>
>> >>
>> >> [[alternative HTML version deleted]]
>> >>
>> >> _______________________________________________
>> >> R-sig-Geo mailing list
>> >> [email protected]
>> >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> >>
>> > [[alternative HTML version deleted]]
>> >
>> > _______________________________________________
>> > R-sig-Geo mailing list
>> > [email protected]
>> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [email protected]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
> _______________________________________________
> R-sig-Geo mailing list
> [email protected]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[[alternative HTML version deleted]]
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
