如何根据实际经纬度绘制地图？

Question

我有一个按经度和纬度划分的区域。我知道1纬度的距离是111,320m，1度经度的距离随纬度的不同而变化。我的网格尺寸是25m×25m，图中所示的网格（总共600个网格）太大了（见下图），不太现实。

我的代码如下：

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"))

我希望网格之间保持适当的间隙，图像的横纵坐标的长度应该更符合实际比例（目前2个纬度对应5个经度，不太实用）。

Answer 1

library(tidyverse)
library(readr)
library(sf)
library(ggspatial)
library(prettymapr)

ggplot_build(p1) %>% 
 .$data %>% .[[1]] %>% 
  add_row(value = color_breaks, group = 0) %>% 
  arrange(desc(group)) %>% 
  fill(xmin:ymax, .direction = "up") %>% 
  mutate(X_1 = xmin, Y_1 = ymin, X_2 = xmin, Y_2 = ymax, 
         X_3 = xmax, Y_3 = ymax, X_4 = xmax, Y_4 = ymin) %>% 
  select(fill, xbin, ybin, value, X_1:Y_4) %>% 
  mutate(poly = row_number()) %>% 
  pivot_longer(cols = X_1:Y_4, values_to = "coord", 
               names_sep = "_", names_to = c("var", "grp")) %>% 
  pivot_wider(names_from = var, values_from = coord) %>% 
  st_as_sf(coords=c("X","Y")) %>% 
  group_by(poly, fill, xbin, ybin, value) %>% 
  summarise() %>% 
  st_cast("POLYGON") %>% 
  st_convex_hull() %>%
  st_set_crs(4326) -> df_plot

ggplot(df_plot) +
  annotation_map_tile(alpha = 0.5) +
  geom_sf(aes(fill = 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
  )),  color = NA) +
  scale_fill_manual(
    values = colors,
    drop = FALSE,
    guide = guide_legend(reverse = TRUE)
  ) +
  labs(
    title = "Data1",
    subtitle = "10 years",
    x = "Longitude (°)",
    y = "Latitude (°)",
    fill = expression("Amount")
  ) +
  theme_minimal() +
  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")
  )

^{创建于 2024-09-03，使用 reprex v2.0.2}