如何基于实际经纬度绘制符合真实比例的网格地图？

我有一个按经纬度划分的区域，已知1度纬度对应111320米，1度经度的距离随纬度变化。目标网格大小为25m×25m，但当前绘制的600个网格过大，不符合实际情况。

我希望保持网格间的合适间距，同时让图像的纵横坐标长度更符合实际比例（当前2度纬度对应5度经度，不符合实际）。

现有代码 #

library(ggplot2)
library(dplyr)
library(readr)
color_breaks <- c(0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0)
colors <- c("#E31A1C", "#FF7F00", "#FDBF6F", "#E9E4A6", "#A4D4A9", "#B2DF8A", "#33A02C", "#1F78B4")

p1 <- ggplot(merged_data, aes(x = longitude, y = latitude, z = median_mafruit)) +
  labs(
    title = "Data1",
    subtitle = "10 years",
    x = "Longitude (°)",
    y = "Latitude (°)",
    fill = expression("Amount")
  ) +
  theme_minimal()

p1<-p1 +
  stat_summary_2d(
    aes(fill = after_stat(
      cut(
        value,
        breaks = c(0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0),
        include.lowest = TRUE,
        right = FALSE
      )
    )),
    bins = 30,
    show.legend = TRUE
  ) +
  scale_fill_manual(
    values = colors,
    drop = FALSE,
    guide = guide_legend(reverse = TRUE)
  ) +
  theme(
    plot.title = element_text(size = 14, face = "bold"),
    plot.subtitle = element_text(size = 12, face = "italic"),
    legend.title = element_text(size = 10),
    legend.text = element_text(size = 8),
    legend.key.size = unit(0.5, "cm"),
    legend.spacing.x = unit(0.5, "cm"),
    legend.spacing.y = unit(0.5, "cm")
  )
p1

部分数据集 #

structure(list(Order = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 
13, 14, 15, 16, 17, 18, 19, 20), latitude = c(43.2143, 43.3697, 
43.3909, 43.3926, 43.3961, 43.3978, 43.066, 43.2215, 43.368, 
43.435, 43.4434, 43.4623, 43.4895, 43.4856, 43.3738, 43.4761, 
43.5102, 43.5118, 43.5062, 43.4933), longitude = c(-4.479, -4.4804, 
-4.4606, -4.4484, -4.4241, -4.412, -4.1046, -4.1049, -4.1031, 
-4.0818, -4.021, -3.9502, -3.8184, -3.7798, -3.7279, -3.7147, 
-3.5971, -3.5849, -3.5585, -3.5177), median_mafruit = c(2.73, 
1.095, 1.115, 2.73, 0.527, 0.527, 0.962, 1.039, 1.039, 2.73, 
2.73, 2.73, 2.73, 2.73, 0.544, 2.73, 2.73, 2.73, 0.478, 2.73)), row.names = c(NA, 
-20L), spec = structure(list(cols = list(Order = structure(list(), class = c("collector_double", 
"collector")), latitude = structure(list(), class = c("collector_double", 
"collector")), longitude = structure(list(), class = c("collector_double", 
"collector")), median_mafruit = structure(list(), class = c("collector_double", 
"collector"))), default = structure(list(), class = c("collector_guess", 
"collector")), delim = ","), class = "col_spec"), class = c("spec_tbl_df", 
"tbl_df", "tbl", "data.frame"))

解决方案 #

1. 核心思路

• 放弃固定bins=30，根据25m目标网格换算出对应的经纬度间隔，动态计算分箱数量
• 使用coord_sf()保证地理坐标的真实比例，避免图形拉伸

2. 修改后的完整代码

library(ggplot2)
library(dplyr)
library(readr)

color_breaks <- c(0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0)
colors <- c("#E31A1C", "#FF7F00", "#FDBF6F", "#E9E4A6", "#A4D4A9", "#B2DF8A", "#33A02C", "#1F78B4")

# 计算25m网格对应的经纬度间隔
center_lat <- mean(merged_data$latitude)
# 1度经度对应的米数：111320 * cos(纬度弧度)
lon_per_meter <- 1 / (111320 * cos(center_lat * pi / 180))
# 1度纬度对应的米数固定为1/111320
lat_per_meter <- 1 / 111320
# 25m对应的经纬度间隔
target_lon_bin <- 25 * lon_per_meter
target_lat_bin <- 25 * lat_per_meter

# 计算覆盖数据范围所需的分箱数
x_bins <- diff(range(merged_data$longitude)) / target_lon_bin
y_bins <- diff(range(merged_data$latitude)) / target_lat_bin

p1 <- ggplot(merged_data, aes(x = longitude, y = latitude, z = median_mafruit)) +
  labs(
    title = "Data1",
    subtitle = "10 years",
    x = "Longitude (°)",
    y = "Latitude (°)",
    fill = expression("Amount")
  ) +
  theme_minimal() +
  # 使用动态计算的分箱数
  stat_summary_2d(
    aes(fill = after_stat(
      cut(
        value,
        breaks = color_breaks,
        include.lowest = TRUE,
        right = FALSE
      )
    )),
    bins = c(x_bins, y_bins),
    show.legend = TRUE
  ) +
  scale_fill_manual(
    values = colors,
    drop = FALSE,
    guide = guide_legend(reverse = TRUE)
  ) +
  # 强制使用地理坐标比例
  coord_sf(crs = 4326) +
  theme(
    plot.title = element_text(size = 14, face = "bold"),
    plot.subtitle = element_text(size = 12, face = "italic"),
    legend.title = element_text(size = 10),
    legend.text = element_text(size = 8),
    legend.key.size = unit(0.5, "cm"),
    legend.spacing.x = unit(0.5, "cm"),
    legend.spacing.y = unit(0.5, "cm")
  )
p1

关键说明

• 网格精度控制：通过中心纬度计算经度的米数换算系数，确保每个网格的实际尺寸接近25m×25m；若数据覆盖纬度范围大，可改为按每个点的纬度单独计算间隔
• 比例修正：coord_sf(crs=4326)会自动按照WGS84坐标系的实际距离比例渲染图形，彻底解决经纬度长度比例不符的问题

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