如何对SpatialPointsFrame与SpatialPolygonsFrame做空间连接并绘制点密度图？

针对你拥有SpatialPointsDataFrame和SpatialPolygonsDataFrame的场景，用sf库可以快速完成空间连接与可视化，步骤如下：

1. 加载依赖库 #

首先加载所需的R包：

library(sf)
library(ggplot2)
library(sp)  # 用于处理原Spatial对象

2. 将Spatial对象转换为sf格式 #

sf库的操作基于其原生格式，先把原sp对象转换为sf对象：

# 转换点数据
points_sf <- st_as_sf(your_spatial_points)  # 替换为你的SpatialPointsDataFrame对象名
# 转换多边形数据
polygons_sf <- st_as_sf(your_spatial_polygons)  # 替换为你的SpatialPolygonsDataFrame对象名

3. 空间连接并统计点数 #

通过空间连接匹配点与所属多边形，再统计每个多边形内的点数：

# 执行点在多边形内的空间连接
point_poly_match <- st_join(points_sf, polygons_sf, join = st_within)

# 统计每个多边形的点数（需替换为你的多边形唯一ID列名，比如"poly_id"）
poly_point_count <- polygons_sf %>%
  left_join(
    point_poly_match %>%
      st_drop_geometry() %>%  # 移除几何信息以高效统计
      count(your_poly_id_column),  # 替换为多边形数据中的唯一标识列
    by = "your_poly_id_column"
  ) %>%
  mutate(n = ifelse(is.na(n), 0, n))  # 无点的多边形点数设为0

4. 绘制可视化图 #

用ggplot2结合sf绘制带颜色映射的多边形分布图：

ggplot(poly_point_count) +
  geom_sf(aes(fill = n), color = "white", linewidth = 0.3) +  # 白色边框区分多边形
  scale_fill_viridis_c(option = "plasma", name = "点数") +  # 色盲友好的配色方案
  theme_minimal() +
  labs(title = "多边形内点数分布", x = "", y = "")

注意事项 #

• 确保点和多边形数据的CRS（坐标参考系统）一致，可通过st_crs(points_sf)和st_crs(polygons_sf)检查，不一致时用st_transform(polygons_sf, st_crs(points_sf))转换。
• 如果你的多边形没有唯一ID列，可先用polygons_sf$poly_id <- seq(nrow(polygons_sf))生成。

内容的提问来源于stack exchange，提问作者Nitzan Sagie