Площадь графика с использованием geom_sf в ggplot

Поскольку координаты ne_countries() определяются географической системой координат (GCS), это затрудняет создание квадратного графика. Чтобы упростить создание квадратного графика, я рекомендую спроецировать ваши данные в систему координат проекции (PCS), например WGS84 Pseudo Mercator/EPSG:3857 или соответствующую зону UTM.

Рабочий процесс:

1. Данные проекта с использованием sf::st_transform()
2. Точки возврата ограничивающей рамки ваших данных
3. Определите, какая широта или долгота больше, и верните половину разницы между ними, например, переменная diff в этом примере.
4. Используйте diff, чтобы увеличить размер меньшего измерения на вашем графике, например. долгота в этом примере

Обратите внимание: в зависимости от того, где находятся ваши данные, разница между исходной версией WGS84/EPSG:3857 и вашей версией PCS может различаться. Однако существует ряд доступных PCS, которые уменьшат любую разницу относительно того, где находятся ваши данные.

library(rnaturalearth)
library(dplyr)
library(ggplot2)
library(patchwork)

# Load ne_countries()
world <- ne_countries(scale = "medium", returnclass = "sf")

# Create projected and unprojected sf objects
italy3857 <- filter(world, sovereignt == "Italy") |> st_transform(3857)
italy4326 <- filter(world, sovereignt == "Italy")

# Function to get polygon from boundary box
bbox_polygon <- function(x) {
  bb <- sf::st_bbox(x)
  
  p <- matrix(
    c(bb["xmin"], bb["ymin"], 
      bb["xmin"], bb["ymax"],
      bb["xmax"], bb["ymax"], 
      bb["xmax"], bb["ymin"], 
      bb["xmin"], bb["ymin"]),
    ncol = 2, byrow = T
  )
  
  sf::st_polygon(list(p))
}

# Create points of bbox
sf_p <- st_sfc(bbox_polygon(italy3857)) |>
  st_as_sf(crs = st_crs(italy3857)) |>
  st_cast("POINT")

sf_p
# Simple feature collection with 5 features and 0 fields
# Geometry type: POINT
# Dimension:     XY
# Bounding box:  xmin: 737796 ymin: 4395685 xmax: 2057834 ymax: 5955490
# Projected CRS: WGS 84 / Pseudo-Mercator
#                         x
# 1  POINT (737796 4395685)
# 2  POINT (737796 5955490)
# 3 POINT (2057834 5955490)
# 4 POINT (2057834 4395685)
# 5  POINT (737796 4395685)

# Get longitude distance
x <- as.integer(st_distance(sf_p[1,], sf_p[4,]))
x
# [1] 1320038

# Get latitude distance
y <- as.integer(st_distance(sf_p[1,], sf_p[2,]))
y
# [1] 1559805

# Subtract max from min, divide by 2
diff <- (y - x) / 2

# Plot unprojected sf (for comparison)
p1 <- ggplot() +
  geom_sf(data = italy4326) +
  labs(title = "WGS84/EPSG:4326")

p2 <- ggplot() +
  geom_sf(data = italy3857) +
  coord_sf(xlim = c(st_bbox(italy3857)[[1]] - diff, 
                    st_bbox(italy3857)[[3]] + diff)) +
  labs(title = "WGS84 Pseudo Mercator/EPSG:3857")

p1 + p2