R地图记

栏目: R语言 · 发布时间: 4年前

内容简介:可灵活调整图形的任意组成成分,同时可在图形上添加2个或多个维度的数据;绘制基础地图方法,仍然只能绘制一维的数据。同时绘制的地图依赖google地图;国人开发的基于百度地图Echart。优点,绘制地图方便快捷,省市级地区的二级地图非常精准,并可绘制炫酷的迁徙图和热图,推荐学习网址:

在重测序文章中经常见到用地图来描述测序样品分布,在R中可轻松复现。

R地图记

1. 不同方法比较

ggplot2

可灵活调整图形的任意组成成分,同时可在图形上添加2个或多个维度的数据;

maps

简单易操作,但原先中国的基础地图中,没有将四川和重庆区分开,现在虽然已经区分,但每个省份轮廓看起来还是与地图略有区别( 国家基础地理信息中心 );

R地图记

googleVis

绘制基础地图方法,仍然只能绘制一维的数据。同时绘制的地图依赖google地图;

REmap

国人开发的基于百度地图Echart。优点,绘制地图方便快捷,省市级地区的二级地图非常精准,并可绘制炫酷的迁徙图和热图,推荐学习网址: http://lchiffon.github.io/REmap/ ;缺点,同googleVis一样,只能绘制一维的数据,同时地图上只能显示中文地名。

2. 地图数据下载

Download GADM data

R地图记 但是从GDM网站下载的中国地图没有台湾,果断差评。

GIS数据

http://cos.name/wp-content/uploads/2009/07/chinaprovinceborderdata_tar_gz.zip

主要是下载三个中国行政区地图数据信息文件: bou2_4p.dbf,bou2_4p.shp和bou2_4p.shx;

使用中如果出现中文省份名称乱码,设置 Sys.setlocale("LC_ALL", "chinese") 即可。

中国行政区地图数据信息数据中包含了925条记录,每条记录中都含有

面积(AREA)

周长(PERIMETER)

各种编号,ADCODE99 是国家基础地理信息中心定义的区域代码,共有 6 位数字,由省、地市、县各两位代码组成。

中文名(NAME)等字段,其中中文名(NAME)字段是以GBK编码的。可利用iconv 格式转换函数来转换各省名称 table(iconv(map$NAME, from = "GBK"))

解压后三个文件放到相同目录下;虽然只读取.shp 文件,.shx 和 .dbf文件也必须在同一个文件目录下才能读取成功。

3. 地图绘制

1. Preparation

setwd("F:/Rwork/china_map")
library(maptools)
library(dplyr)
library(ggplot2)
library(RColorBrewer)
Sys.setlocale("LC_ALL", "chinese")

2. Map Data

Download GIS数据: http://cos.name/wp-content/uploads/2009/07/chinaprovinceborderdata_tar_gz.zip

解压后三个文件放到当前目录( getwd() )下;

虽然只读取.shp 文件,.shx 和 .dbf文件也必须在同一个文件目录下才能读取成功。

map_data <- readShapePoly("bou2_4p.shp")
names(map_data)
map_data@data$ID <- row.names(map_data@data)
# 去掉包含NA的数据
map_data@data <- na.omit(map_data@data)
nrow(map_data@data)
# 可选,按照省份面积(AREA)筛选,主要为去掉南沙群岛和围绕南海的许多小岛
Fmap_data <- subset(map_data, AREA > 0.005)
nrow(Fmap_data@data)

添加省会拼音

Create a data.frame called cnmapdf which contains id, prov_en and prov_cn and key map plotting info;

prov_cn <- unique(map_data$NAME)
prov_en <- c("Heilongjiang", "Inner Mongolia", "Xinjiang", "Jilin",
             "Liaoning", "Gansu", "Hebei", "Beijing", "Shanxi",
             "Tianjin", "Shaanxi", "Ningxia", "Qinghai", "Shandong",
             "Tibet", "Henan", "Jiangsu", "Anhui", "Sichuan", "Hubei",
             "Chongqing", "Shanghai", "Zhejiang", "Hunan", "Jiangxi",
             "Yunnan", "Guizhou", "Fujian", "Guangxi", "Taiwan", 
             "Guangdong", "Hong Kong", "Hainan")

prov <- data.frame(prov_cn, prov_en)
id_prov <- map_data@data %>%
  mutate(prov_en = sapply(NAME, function(x) prov$prov_en[which(prov_cn == x)])) %>%
  mutate(prov_cn = as.character(NAME),prov_en = as.character(prov_en)) %>%
  select(id = ID, prov_cn, prov_en)

cnmapdf <- plyr::join(fortify(map_data), id_prov, by = "id")
head(cnmapdf)

添加省会城市坐标

cap_coord <- c(
  "Beijing", "北京", "Beijing", 116.4666667, 39.9,
  "Shanghai", "上海", "Shanghai", 121.4833333, 31.23333333,
  "Tianjin", "天津", "Tianjin", 117.1833333, 39.15,
  "Chongqing", "重庆", "Chongqing", 106.5333333, 29.53333333,
  "Harbin", "哈尔滨", "Heilongjiang", 126.6833333, 45.75,
  "Changchun", "长春", "Jilin", 125.3166667, 43.86666667,
  "Shenyang", "沈阳", "Liaoning", 123.4, 41.83333333,
  "Hohhot", "呼和浩特", "Inner Mongolia", 111.8, 40.81666667,
  "Shijiazhuang", "石家庄", "Hebei", 114.4666667, 38.03333333,
  "Taiyuan", "太原", "Shanxi", 112.5666667, 37.86666667,
  "Jinan", "济南","Shandong", 117, 36.63333333,
  "Zhengzhou", "郑州", "Henan", 113.7, 34.8, 
  "Xi'an", "西安", "Shaanxi", 108.9, 34.26666667,
  "Lanzhou", "兰州", "Gansu", 103.8166667, 36.05,
  "Yinchuan", "银川", "Ningxia", 106.2666667, 38.33333333,
  "Xining", "西宁", "Qinghai", 101.75, 36.63333333,
  "Urumqi", "乌鲁木齐", "Xinjiang", 87.6, 43.8,
  "Hefei", "合肥", "Anhui", 117.3, 31.85,
  "Nanjing", "南京", "Jiangsu", 118.8333333, 32.03333333,
  "Hangzhou", "杭州", "Zhejiang", 120.15, 30.23333333,
  "Changsha", "长沙", "Hunan", 113, 28.18333333,
  "Nanchang", "南昌", "Jiangxi", 115.8666667, 28.68333333,
  "Wuhan", "武汉", "Hubei", 114.35, 30.61666667,
  "Chengdu", "成都", "Sichuan", 104.0833333, 30.65,
  "Guiyang", "贵阳", "Guizhou", 106.7, 26.58333333,
  "Fuzhou", "福州", "Fujian", 119.3, 26.08333333,
  "Taibei", "台北", "Taiwan", 121.5166667, 25.05,
  "Guangzhou", "广州", "Guangdong", 113.25, 23.13333333,
  "Haikou", "海口", "Hainan", 110.3333333, 20.03333333,
  "Nanning", "南宁", "Guangxi", 108.3333333, 22.8,
  "Kunming", "昆明", "Yunnan", 102.6833333, 25,
  "Lhasa", "拉萨", "Tibet", 91.16666667, 29.66666667,
  "Hong Kong", "香港", "Hong Kong", 114.1666667, 22.3,
  "Macau", "澳门", "Macau", 113.5, 22.2)

cap_coord <- as.data.frame(matrix(cap_coord, nrow = 34, byrow = TRUE))
names(cap_coord) <- c("city_en", "city_cn", "prov_en", "long", "lat")
cap_coord <- cap_coord %>%
  mutate(prov_en = as.vector(prov_en),
         city_en = as.vector(city_en),
         city_cn = as.vector(city_cn),
         cap_long = as.double(as.vector(long)),
         cap_lat = as.double(as.vector(lat))) %>%
  select(prov_en, city_en, city_cn, cap_long, cap_lat)

head(cap_coord)
cnmapdf <- plyr::join(cnmapdf, cap_coord, by = "prov_en", type = "full")

3. 开始绘制地图

选择一个省画图

默认绘制的地图的形状有些扁平。这是因为,在绘图的过程中,默认把经度和纬度作为普通数据,均匀平等对待,绘制在笛卡尔坐标系上造成的。其实,地球的球面图形如何映射到平面图上,在地理学上是有一系列不同的专业算法的。地图不应该画在普通的笛卡尔坐标系上,而是要画在地理学专业的坐标系上。在这一点上,R 的 ggplot2 包提供了专门的coord_map()函数。

shanghai <- cnmapdf[cnmapdf$prov_en == "Shanghai",]
shanghai %>%
  ggplot(aes(x = long, y = lat, group = group, fill=factor(prov_en))) +
  geom_polygon( color = "grey") +
  coord_map() +
  ggtitle("上海直辖市") +
  xlab("经度") +
  ylab("维度") +
  scale_fill_brewer(palette="Paired")

R地图记

画多个省

map1 <- cnmapdf %>%
  filter(prov_en %in% c("Jiangsu", "Zhejiang", "Shanghai"))  %>%
  ggplot() +
  geom_polygon(aes(x = long, y = lat, group = group, fill = prov_cn), color = "grey")

coord_delta_cap <- subset(cap_coord, prov_en %in% c("Zhejiang", "Shanghai", "Jiangsu"))
map1 +
  geom_point(data = coord_delta_cap, aes(x = cap_long, y = cap_lat)) +
  geom_text(data = coord_delta_cap, aes(cap_long, cap_lat - .25, label = city_cn)) +
  coord_map() +
  ggtitle("长江三角洲") +
  xlab("经度") +
  ylab("维度") +
  scale_fill_brewer(palette="Set2")

R地图记

全国地图

map0 <- cnmapdf %>%
  filter(prov_en %in% unique(cnmapdf$prov_en))  %>%
  ggplot() +
  geom_polygon(aes(x = long, y = lat, group = group, fill = "white"), color = "grey") +
  scale_fill_identity()
coord_delta_cap <- subset(cap_coord, prov_en %in% unique(cnmapdf$prov_en))
# 解决重叠地名
library(ggrepel)
spec.city <- c("香港","澳门")
cap_map_data01 <- coord_delta_cap[coord_delta_cap$city_cn %in% spec.city,]
cap_map_data02 <- coord_delta_cap[!coord_delta_cap$city_cn %in% spec.city,]
cnmap <-  map0 + geom_point(data=cap_map_data02,aes(x=cap_long, y= cap_lat),shape=1,colour="white") +
          geom_text(data=cap_map_data02,aes(x=cap_long, y= cap_lat,label=city_cn)) +
          geom_text_repel(data=cap_map_data01,aes(x=cap_long, y= cap_lat,label=city_cn)) +
          coord_map() +
          theme_void() + 
          theme(legend.position = "none") +
          scale_fill_identity()
cnmap

R地图记

map1 <- cnmapdf %>%
  filter(prov_en %in% unique(cnmapdf$prov_en))  %>%
  ggplot() +
  geom_polygon(aes(x = long, y = lat, group = group, fill = prov_cn), color = "grey")

coord_delta_cap <- subset(cap_coord, prov_en %in% unique(cnmapdf$prov_en))

nb.cols <- length(unique(coord_delta_cap$prov_en))
mycolors <- colorRampPalette(brewer.pal(8, "Set2"))(nb.cols)

map1 +
  geom_point(data = coord_delta_cap, aes(x = cap_long, y = cap_lat)) +
  geom_text(data = coord_delta_cap, aes(cap_long, cap_lat - .25, label = city_cn)) +
  coord_map() +
  ggtitle("中国") +
  xlab("经度") +
  ylab("维度") +
  theme_void() + 
  theme(legend.position = "none") +
  scale_fill_manual(values = mycolors)

R地图记

解决重叠地名

library(ggrepel)
spec.city <- c("香港","澳门")
cap_map_data1 <- coord_delta_cap[coord_delta_cap$city_cn %in% spec.city,]
cap_map_data2 <- coord_delta_cap[!coord_delta_cap$city_cn %in% spec.city,]
map1 + geom_point(data=cap_map_data2,aes(x=cap_long, y= cap_lat),shape=1,colour="white") +
  geom_text(data=cap_map_data2,aes(x=cap_long, y= cap_lat,label=city_cn)) +
  geom_text_repel(data=cap_map_data1,aes(x=cap_long, y= cap_lat,label=city_cn)) +
  coord_map() +
  theme_void() + 
  theme(legend.position = "none") +
  scale_fill_manual(values = mycolors)

R地图记

颜色标注全国地图某几个省

https://cosx.org/2009/07/drawing-china-map-using-r/

par(mar=rep(0,4))
library(maps)
library(mapdata)
getColor = function(mapdata, provname, provcol, othercol){
	f = function(x, y) ifelse(x %in% y, which(y == x), 0)
	colIndex = sapply(mapdata@data$NAME, f, provname)
	fg = c(othercol, provcol)[colIndex + 1]
	return(fg)
}
provname = c("北京市", "天津市", "上海市", "重庆市")
provcol = c("red", "green", "yellow", "purple")
plot(map_data, col = getColor(map_data, provname, provcol, "white"))
points(cap_coord$cap_long, cap_coord$cap_lat, pch = 19, col = rgb(0, 0, 0, 0.5))
text(cap_coord$cap_long, cap_coord$cap_lat, cap_coord[, 3], cex = 0.9, col = rgb(0,0, 0, 0.7), 
     pos = c(2, 4, 4, 4, 3, 4, 2, 3, 4, 2, 4, 2, 2, 4, 3, 2, 1, 3, 1, 1, 2, 3, 2, 2, 1, 2, 4, 3, 1, 2, 2, 4, 4, 2))
axis(1, lwd = 0); axis(2, lwd = 0); axis(3, lwd = 0); axis(4, lwd = 0)
as.character(na.omit(unique(map_data@data$NAME)))

R地图记

颜色标注全国地图某几个省 (推荐)

provname = c("北京市", "天津市", "上海市", "重庆市")
provcol = c("red", "green", "yellow", "purple")
getColors = function(mapdata, provname, provcol, othercol){
	f = function(x, y) ifelse(x %in% y, which(y == x), 0)
	colIndex = sapply(mapdata$prov_cn, f, provname)
	fg = c(othercol, provcol)[colIndex + 1]
	return(fg)
}
mc=getColors(cnmapdf, provname, provcol, "white")
map2 <- cnmapdf %>%
  filter(prov_en %in% unique(cnmapdf$prov_en))  %>%
  ggplot() +
  geom_polygon(aes(x = long, y = lat, group = group, fill = mc), color = "grey")

coord_delta_cap <- subset(cap_coord, prov_en %in% unique(cnmapdf$prov_en))
# 解决重叠地名
library(ggrepel)
spec.city <- c("香港","澳门")
cap_map_data21 <- coord_delta_cap[coord_delta_cap$city_cn %in% spec.city,]
cap_map_data22 <- coord_delta_cap[!coord_delta_cap$city_cn %in% spec.city,]
map2 + geom_point(data=cap_map_data22,aes(x=cap_long, y= cap_lat),shape=1,colour="white") +
  geom_text(data=cap_map_data22,aes(x=cap_long, y= cap_lat,label=city_cn)) +
  geom_text_repel(data=cap_map_data21,aes(x=cap_long, y= cap_lat,label=city_cn)) +
  coord_map() +
  theme_void() + 
  theme(legend.position = "none") +
  scale_fill_identity()

R地图记

4. 地图添加数据

实例数据下载: 中华人民共和国国家统计局

Heatmap

democn <- read.csv("China_pop.csv", stringsAsFactors = F, check.names=FALSE)
library(tidyr)
library(reshape2)
democndf <- melt(democn,variable.name ="year", value.name = "population")
head(spread(democndf, year, population))

map2df <- cnmapdf %>% 
  plyr::join(subset(democndf, year == "2018年"), by = "prov_cn") %>%
  mutate(population = as.numeric(population))
  
map2df %>%
  ggplot() +
  geom_polygon(aes(x = long, y = lat, group = group, fill = population), color = "grey") +
  geom_point(data=cap_map_data22,aes(x=cap_long, y= cap_lat),shape=1,colour="white") +
  geom_text(data=cap_map_data22,aes(x=cap_long, y= cap_lat,label=city_cn),size=2) +
  geom_text_repel(data=cap_map_data21,aes(x=cap_long, y= cap_lat,label=city_cn),size=2) +
  scale_fill_gradient(low = "red", high = "yellow") +
  theme_void()

R地图记 多个图

map3df <- cnmapdf %>% 
  plyr::join(democndf, by = "prov_cn") %>%
  mutate(population = as.numeric(population)) %>%
  na.omit()
map3df %>%
  ggplot(aes(x = long, y = lat, group = group, fill = population)) +
  geom_polygon(color = "grey", lwd = .1) +
  coord_equal() +
  facet_wrap(~year)

R地图记

Bubbles

map1 + 
  geom_point(data = map2df, aes(cap_long, cap_lat, size = population), shape = 21, fill="#9070c7",colour="grey", alpha = .5) +
  scale_size_area(max_size=5) +
  geom_text(data=cap_map_data2,aes(x=cap_long, y= cap_lat,label=city_cn),size=2,vjust=0,nudge_y=0.5) +
  geom_text_repel(data=cap_map_data1,aes(x=cap_long, y= cap_lat,label=city_cn),size=2,vjust=0,nudge_y=0.5) +
  coord_map() +
  theme_void() + 
  theme(legend.position = "none") +
  scale_fill_manual(values = mycolors)

R地图记

Bar

map1 + 
  geom_errorbar(data=map2df,aes(x=cap_long, ymin=cap_lat, ymax=cap_lat + population/3000 ),
                colour="blue",size=2, width=0,alpha=0.5) +
  geom_text(data=cap_map_data2,aes(x=cap_long, y= cap_lat,label=city_cn),size=2,vjust=0,nudge_y=0.5) +
  geom_text_repel(data=cap_map_data1,aes(x=cap_long, y= cap_lat,label=city_cn),size=2,vjust=0,nudge_y=0.5) +
  coord_map() +
  theme_void() + 
  theme(legend.position = "none") +
  scale_fill_manual(values = mycolors)

R地图记

5. 世界地图

library(rworldmap)
met <- as.data.frame(read.csv("MetObjects_5k-sample.csv"))
countries.met <- as.data.frame(table(met$Country))
head(countries.met)
colnames(countries.met) <- c("country", "value")
matched <- joinCountryData2Map(countries.met, joinCode="NAME", nameJoinColumn="country")
mapCountryData(matched, nameColumnToPlot="value", mapTitle="Met Collection Country Sample", catMethod = "pretty", colourPalette = "heat",oceanCol="aliceblue")

R地图记

仅显示某一区域

mapCountryData(matched, nameColumnToPlot="value", mapTitle="Met Collection in Eurasia", mapRegion="Eurasia", colourPalette="heat", catMethod="pretty", oceanCol="aliceblue")

R地图记

library(ggplot2)
library(dplyr)

WorldData <- map_data('world') %>% filter(region != "Antarctica") %>% fortify

df <- data.frame(region=c('Hungary','Lithuania','Argentina'), 
                 value=c(4,10,11), 
                 stringsAsFactors=FALSE)

p <- ggplot() +
    geom_map(data = WorldData, map = WorldData,
                  aes(x =long , y = lat, group = group, map_id=region),
                  fill = "white", colour = "#7f7f7f", size=0.5) + 
    geom_map(data = df, map=WorldData, aes(fill=value, map_id=region),colour="#7f7f7f", size=0.5) +
    coord_map("rectangular", lat0=0, xlim=c(-180,180), ylim=c(-60, 90)) +
    scale_fill_continuous(low="thistle2", high="darkred", guide="colorbar") +
    scale_y_continuous(breaks=c()) +
    scale_x_continuous(breaks=c()) +
    labs(fill="legend", title="Title", x="", y="") +
    theme(text = element_text( color = "#FFFFFF")
        ,panel.background = element_rect(fill = "aliceblue")
        ,plot.background = element_rect(fill = "aliceblue")
        ,panel.grid = element_blank()
        ,plot.title = element_text(size = 30)
        ,plot.subtitle = element_text(size = 10)
        ,axis.text = element_blank()
        ,axis.title = element_blank()
        ,axis.ticks = element_blank()
        ,legend.position = "right"
        )
    #theme_bw()
p

R地图记

6. Info

## R version 3.4.2 (2017-09-28)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 17134)
## 
## Matrix products: default
## 
## locale:
## [1] LC_COLLATE=Chinese (Simplified)_China.936 
## [2] LC_CTYPE=Chinese (Simplified)_China.936   
## [3] LC_MONETARY=Chinese (Simplified)_China.936
## [4] LC_NUMERIC=C                              
## [5] LC_TIME=Chinese (Simplified)_China.936    
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] rworldmap_1.3-6    reshape2_1.4.3     tidyr_0.8.3       
##  [4] mapdata_2.3.0      maps_3.3.0         ggrepel_0.8.0     
##  [7] RColorBrewer_1.1-2 ggplot2_3.1.0      dplyr_0.8.0.1     
## [10] maptools_0.9-5     sp_1.3-1          
## 
## loaded via a namespace (and not attached):
##  [1] Rcpp_1.0.0       pillar_1.3.1     compiler_3.4.2   plyr_1.8.4      
##  [5] tools_3.4.2      dotCall64_1.0-0  digest_0.6.18    evaluate_0.13   
##  [9] tibble_2.0.1     gtable_0.2.0     lattice_0.20-38  pkgconfig_2.0.2 
## [13] rlang_0.3.1      mapproj_1.2.6    yaml_2.2.0       spam_2.2-2      
## [17] xfun_0.5         withr_2.1.2      stringr_1.4.0    knitr_1.21      
## [21] fields_9.8-3     grid_3.4.2       tidyselect_0.2.5 glue_1.3.0      
## [25] R6_2.4.0         foreign_0.8-71   rmarkdown_1.11   purrr_0.3.1     
## [29] searcher_0.0.3   magrittr_1.5     scales_1.0.0     htmltools_0.3.6 
## [33] assertthat_0.2.0 colorspace_1.4-0 labeling_0.3     stringi_1.3.1   
## [37] lazyeval_0.2.1   munsell_0.5.0    crayon_1.3.4

4. 参考

R Visual. - China Map Part II

https://www.datanovia.com/en/blog/ggplot-colors-best-tricks-you-will-love/

https://www.datanovia.com/en/blog/top-r-color-palettes-to-know-for-great-data-visualization/


以上就是本文的全部内容,希望本文的内容对大家的学习或者工作能带来一定的帮助,也希望大家多多支持 码农网

查看所有标签

猜你喜欢:

本站部分资源来源于网络,本站转载出于传递更多信息之目的,版权归原作者或者来源机构所有,如转载稿涉及版权问题,请联系我们

JSP应用开发技术

JSP应用开发技术

柳永坡 / 人民邮电出版社 / 2005-9 / 52.00元

本书全面系统地介绍了JSP应用开发技术,包括JSP预备知识和环境配置、JSP编程基础、JSP应用开发进阶、在JSP中使用数据库、Servlet技术、标签库和表达式语言、Web编程模式和应用框架等几个方面的内容。本书不但由浅入深地介绍了JSP程序设计的原理、方法和技术,还提供了大量的JSP应用开发实例,给出了相应的实用技巧、操作步骤及优化思路。 本书着重于JSP技术的应用性和可操作性,......一起来看看 《JSP应用开发技术》 这本书的介绍吧!

MD5 加密
MD5 加密

MD5 加密工具

Markdown 在线编辑器
Markdown 在线编辑器

Markdown 在线编辑器

UNIX 时间戳转换
UNIX 时间戳转换

UNIX 时间戳转换