1、圆柱投影地图
圆柱投影是将一个圆柱面包围椭球体，并使之相切或相割，再根据某种条件将椭球面上的经纬网点投影到圆柱面上，然后，沿圆柱面的一条母线切开，将其展成平面而得到的投影。圆柱投影地图往往和现实差别太大，陆地变形非常严重，离赤道越远的形状就越膨胀，特别在南极和北极，比如圆柱投影地图上的格陵兰岛的面积看起来与非洲类似（实际上只有非洲的1/4）。常见的圆柱投影有墨卡托、横轴墨卡托和米勒投影。
## 加载R包
library(rgdal)
library(ggplot2)
library(maps)
library(RColorBrewer)
## 读取maps中世界地图国家边界数据
world_map <- map_data("world")
##绘制世界地图
ggplot(world_map, aes(x = long, y = lat, group = group)) +
geom_polygon(fill = "gray50", colour = "gray50", linetype = "dashed") +
coord_equal() +
theme_bw()
2、伪圆柱投影地图
伪圆柱投影是在圆柱投影基础上，规定纬线为平行直线，而经线则根据某些特定条件而设计成对称于中央经线的各类曲线（多为正弦曲线或椭圆曲线）的投影。常见的伪圆柱投影有埃克特投影、古蒂等面积投影、自然地球 (Natural Earth) 投影、罗宾森 (Robinson) 投影、瓦格纳投影、摩尔维德 (Mollweide) 投影、温克尔投影等。罗宾森投影可能是一种最常用于绘制世界地图的折衷伪圆柱地图投影。“国家地理”将罗宾森投影用于其世界地图约十年，直到 1998 年。1998年后，温克尔三重 (Winkel Tripel) 投影常被美国国家地理学会用于一般世界地图。
# 加载R包
library(sf)
library(rnaturalearth)
library(rnaturalearthdata)
# 读取世界地图
world <- ne_coastline(returnclass = "sf")
# 重置坐标系
eckertIV <- "+proj=eck4 +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs"
# 绘制地图
ggplot(world) +
geom_sf(color = "gray50") +
coord_sf(crs = eckertIV) +
theme_bw()
# 读取地图数据
grat <- readOGR("ne_110m_graticules_all", layer = "ne_110m_graticules_15")
bbox <- readOGR("ne_110m_graticules_all", layer = "ne_110m_wgs84_bounding_box")
# 重置坐标系
grat_robin <- spTransform(grat, CRS("+proj=robin"))
grat_df_robin <- fortify(grat_robin)
bbox_robin <- spTransform(bbox, CRS("+proj=robin"))
bbox_robin_df <- fortify(bbox_robin)
# 绘图——Robinson投影地图
ggplot(bbox_robin_df, aes(long, lat, group = group)) +
geom_polygon(fill = "white") +
geom_polygon(data = wmap_df_robin, aes(long, lat, group = group, fill = hole)) +
geom_path(data = grat_df_robin, aes(long, lat, group = group, fill = NULL), linetype = "dashed", color = "grey50") +
labs(title = "World map") +
coord_equal() +
theme_bw() +
scale_fill_manual(values = c("gray50", "white"), guide = "none")
# 读取地图文件
grat <- readOGR("ne_110m_graticules_all", layer = "ne_110m_graticules_15")
# 重置坐标系
grat_robin <- spTransform(grat, CRS("+proj=robin"))
grat_df_robin <- fortify(grat_robin)
# add country borders
countries <- readOGR("ne_110m_admin_0_countries", layer = "ne_110m_admin_0_countries")
# countries_robin <- spTransform(countries, CRS("+init=ESRI:54030")) ## not run
countries_robin <- spTransform(countries, CRS("+proj=robin"))
countries_robin_df <- fortify(countries_robin)
# 读取数据
places_df <- read.csv("sample location.csv") ## Contains longitude, latitude information and sample size
# 选取对应的经纬度信息和分组信息
places_robin_df <- project(cbind(places_df$longitude, places_df$latitude), proj = "+proj=robin")
places_robin_df <- as.data.frame(places_robin_df)
names(places_robin_df) <- c("LONGITUDE", "LATITUDE")
# 分组信息
places_robin_df$samples <- places_df$scalerank
ggplot(bbox_robin_df, aes(long, lat, group=group)) +
geom_polygon(fill = "white") +
geom_polygon(data = countries_robin_df, aes(long, lat, group = group, fill = hole)) +
geom_point(data = places_robin_df, aes(LONGITUDE, LATITUDE, group = NULL, fill = NULL, size = samples), color = "#32caf6", alpha = I(8/10)) +
geom_path(data = countries_robin_df, aes(long, lat, group = group, fill = hole), color = "white", size = 0.3) +
geom_path(data = grat_df_robin, aes(long, lat, group = group, fill = NULL), linetype = "dashed", color = "grey50") +
labs(title = "World map") +
coord_equal() +
theme_bw() +
scale_fill_manual(values = c("black", "white"), guide = "none")+
scale_size_continuous(range = c(1,4), guide = "none")
# 读取数据
countries <- readOGR("ne_110m_admin_0_countries", layer = "ne_110m_admin_0_countries")
bbox <- readOGR("ne_110m_graticules_all", layer = "ne_110m_wgs84_bounding_box")
grat <- readOGR("ne_110m_graticules_all", layer = "ne_110m_graticules_15")
wmap <- readOGR("ne_110m_land", layer = "ne_110m_land")
# 重置坐标系
countries_wintri <- spTransform(countries, CRS("+proj=wintri"))
bbox_wintri <- spTransform(bbox, CRS("+proj=wintri"))
wmap_wintri <- spTransform(wmap, CRS("+proj=wintri"))
grat_wintri <- spTransform(grat, CRS("+proj=wintri"))
## 绘制地图——Winkel三重投影地图
ggplot(bbox_wintri, aes(long,lat, group = group)) +
geom_polygon(fill = "white") +
geom_polygon(data = countries_wintri, aes(long, lat, group = group, fill = hole)) +
geom_path(data = countries_wintri, aes(long, lat, group = group, fill = hole), color = "white", size = 0.3) +
geom_path(data = grat_wintri, aes(long, lat, group = group, fill = NULL), linetype = "dashed", color = "grey50") +
labs(title = "World map") +
coord_equal(ratio = 1) +
theme_bw() +
scale_fill_manual(values = c("black", "white"), guide = "none")
转自：“生态科研笔记”微信公众号
如有侵权，请联系本站删除！
