sf with sf list column

sf_sf_list_column.R

library(mapview)
library(sf)

## create polygon index column
franconia$index = sprintf("%02.f", 1:nrow(franconia))

## join polygons to points
brew_franc = st_join(breweries, franconia[, c("index", "geometry")])

## split on polygon entity
brew_list = split(brew_franc, brew_franc$index)

## include list in polygon data as list column
franc = franconia[franconia$index %in% names(brew_list), ]
franc$breweries = I(brew_list)

## view selected feature
mapview(franc[4, ]) + mapview(franc[4, ]$breweries[[1]])