各纬度海洋面积占比及陆地形态复杂度计算技术问询

一、各纬度海洋占比计算与折线图绘制（R实现） #

核心逻辑 #

将地球按纬度划分为1°间隔的环形带，计算每个带内的海洋面积（带总面积减去陆地面积），再除以带总面积得到占比。必须使用等面积投影保证面积计算精准，这里推荐LAEA（Lambert Azimuthal Equal-Area）投影。

步骤与代码 #

1. 安装并加载依赖包

install.packages(c("rnaturalearth", "sf", "dplyr", "ggplot2"))
library(rnaturalearth)
library(sf)
library(dplyr)
library(ggplot2)

2. 获取世界陆地数据（包含内陆海边界）

# 中等分辨率陆地数据，兼顾精度与计算速度
world_land <- ne_countries(scale = "medium", returnclass = "sf")

3. 生成纬度带多边形并计算占比

# 创建单个纬度带的函数（±0.5°范围，避免极点异常）
create_latitude_band <- function(lat) {
  lat_min <- lat - 0.5
  lat_max <- lat + 0.5
  # 生成环形多边形并转换到等面积投影
  st_polygon(list(rbind(
    c(-180, lat_min), c(180, lat_min), c(180, lat_max), c(-180, lat_max), c(-180, lat_min)
  ))) %>%
    st_sfc(crs = 4326) %>%
    st_transform(crs = "+proj=laea +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs")
}

# 遍历所有纬度（-89°到89°，步长1°）
latitudes <- seq(-89, 89, by = 1)
lat_band_stats <- lapply(latitudes, function(lat) {
  band <- create_latitude_band(lat)
  band_area <- as.numeric(st_area(band))
  # 裁剪陆地到当前纬度带并计算面积
  land_in_band <- st_intersection(world_land, band)
  land_area <- as.numeric(sum(st_area(land_in_band)))
  # 计算海洋占比
  ocean_ratio <- (band_area - land_area) / band_area
  data.frame(latitude = lat, ocean_ratio = ocean_ratio)
}) %>% bind_rows()

4. 绘制折线图

ggplot(lat_band_stats, aes(x = latitude, y = ocean_ratio)) +
  geom_line(color = "#2c7fb8", linewidth = 1) +
  scale_y_continuous(labels = scales::percent_format(accuracy = 1)) +
  labs(x = "纬度", y = "海洋占该纬度带表面积比例", title = "各纬度海洋面积占比") +
  theme_minimal()

二、各纬度陆地形态复杂度量化（R实现） #

核心逻辑 #

在每个纬度绘制一条横跨全球的直线，统计直线与陆地边界的交点数量——交点越多，说明该纬度陆地形态越破碎复杂。同样需要在等面积投影下进行计算，避免交点偏移。

步骤与代码 #

1. 提取陆地边界线

# 将陆地多边形转换为多线段，用于交点计算
world_land_lines <- st_cast(world_land, "MULTILINESTRING")

2. 生成纬度线并统计交点

# 创建单个纬度线的函数
create_latitude_line <- function(lat) {
  st_linestring(rbind(c(-180, lat), c(180, lat))) %>%
    st_sfc(crs = 4326) %>%
    st_transform(crs = "+proj=laea +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs")
}

# 计算每个纬度的交点数（去重避免重复统计角点）
lat_complexity <- lapply(latitudes, function(lat) {
  line <- create_latitude_line(lat)
  intersections <- st_intersection(line, world_land_lines)
  # 提取唯一交点坐标并计数
  unique_points <- unique(st_coordinates(intersections)[,1:2])
  intersection_count <- nrow(unique_points)
  data.frame(latitude = lat, intersection_count = intersection_count)
}) %>% bind_rows()

3. 绘制复杂度折线图

ggplot(lat_complexity, aes(x = latitude, y = intersection_count)) +
  geom_line(color = "#d95f02", linewidth = 1) +
  labs(x = "纬度", y = "陆地边界交点数", title = "各纬度陆地形态复杂度") +
  theme_minimal()

关键注意事项 #

• 投影选择：所有空间计算必须在等面积投影下进行，WGS84等地理坐标系无法精准计算面积和交点。
• 数据精度：若需要更高精度，可将ne_countries的scale参数改为"large"，但计算时间会增加。
• 极点处理：避免直接计算90°N/S，因为极点是单点，无法形成环形带，所以只计算到89°。

备选Python方案 #

使用geopandas、matplotlib实现，核心逻辑与R一致：

1. 用geopandas.read_file加载世界陆地数据
2. 生成纬度带/纬度线，通过geopandas.overlay或intersection进行空间计算
3. 用matplotlib.pyplot绘制折线图

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