library(leaflet)
library(sf)
library(dplyr)
districts <- read_sf(
  system.file("extdata/districts.geojson", package = "geospaar")
)
bb <- unname(st_bbox(districts))
xy <- st_centroid(districts) %>% st_coordinates() %>% 
  bind_cols(name = districts$distName, .)
slist <- list("color" = "white")
label_opts <- labelOptions(noHide = TRUE, style = slist, direction = 'center',
                           textOnly = TRUE, textsize = "5px")

m <- leaflet() %>% 
  addProviderTiles("Esri.WorldImagery") %>% 
  fitBounds(bb[1], bb[2], bb[3], bb[4]) %>% 
  addPolygons(data = districts, fill = FALSE, color = "white", 
              group = "Districts", weight = 1) %>% 
  addLabelOnlyMarkers(xy$X, xy$Y, label = xy$name, group = "Names",
                      labelOptions = label_opts) %>%
  addLayersControl(overlayGroups = c("Districts", "Names"),
                   options = layersControlOptions(collapsed = FALSE,
                                                  autoZIndex = FALSE))

Split-apply-combine

set.seed(10)
dat <- data.frame(v1 = 1:100, v2 = sample(0:10, size = 100, replace = TRUE), 
                  grp = sample(letters[1:3], size = 100, replace = TRUE))

sapply(sort(unique(dat$grp)), function(x) {
  colMeans(dat[dat$grp == x, c("v1", "v2")])
})

dat %>% 
  group_by(grp) %>% 
  summarise(across(c(v1, v2), mean))

Reading, reshaping, analyzing, plotting

dat <- readr::read_csv(
  system.file("extdata/FAOSTAT_maize.csv", package = "geospaar")
)