我找到的解决方案来自于Jared Lander的一个视频，他在视频中讨论了使用h3jsr h3六边形库作为基础六边形形状。数据被转换为空间对象，并与H3地址匹配。任何分析都可以在这个hexbin（平均值/总和）内完成，然后可以在leaflet中绘制。

library(dplyr)
library(h3jsr)
library(leafgl) #for faster plotting of hexbins
library(sf)

n <- 100000
df <- data.frame(lat = rnorm(n, 42.0285, 1),
                 lng = rnorm(n, -93.65, 1),
                 abun = runif(n, 0, 10000))

#dataframe to spatial object
df2 <- st_as_sf(x = df,                         
                coords = c("lng", "lat"),
                crs = "+datum=WGS84",
                na.fail=FALSE)

#Creates the hexbins for the coverage area and then calculate mean abun in each hexbin
# res=4 controls the size of the hexbin; smaller # = larger hexbins.  
hex5 <-df2 |>
  mutate(hex=point_to_cell(geometry, res=4)) |> #geometry column matches to H3 address
  st_drop_geometry()  |> #drop the geometry for faster processing
  group_by(hex) %>% #mean abundance by hex 
  summarise(mean_abun=mean(abun),
            .groups = 'drop')%>% 
  cell_to_polygon(simple=FALSE)#back to spatial 


#mapping the hexbins
leaflet("map") %>%
  addProviderTiles('CartoDB.Voyager')%>%addScaleBar(position='bottomleft') %>%
  addGlPolygons(data=hex5, fillColor="mean_abun",popup = ~mean_abun)