由于近期疫情的播散,而流行病学研究对疫情防控又至关重要,所以,最近涌现了一大批关于疾病流调的文章,这也使得很多研究人员在文章中需要绘制不同地区的地图作为文章中的主图。这些图频繁出现在Lancet或者NEJM、CNS等顶级杂志中,不过对于很多科研人员来说,地图的绘制由于没有现成软件可以直接操作,一直以来都是一大难题。
前面一期的教程中,我们给大家讲解了中国地图的绘制方法【科研猫·绘图】中国地图绘制:R语言版,在我们发布这篇教程后,马上有同学反映想要学习世界地图的绘制,作为实力宠粉的科研猫团队,粉丝的要求岂能不去满足。上期教大家绘制中国地图还不够,因为很多疾病(包括慢性病,不光是传染病),它们的范围是世界性的,所以确实需要学习一下如何绘制世界地图。闲话到此,直奔主题。
LEVEL 1
绘制地图之前,需要获取世界地图的数据,这个数据我们通过R包maptools可以进行读取。不过我们细心的技术在使用这个R包的过程中发现了一个原则性问题,该R包中共有246个国家的名称,其中我国的Taiwan省被列为一个独立的国家。这种事情是坚决不给允许的,于是乎,我们在给大家的文件中删除了Taiwan,并且在代码上进行了调整,使得Taiwan的数据永远与China保持一致,严格贯彻“一个中国”原则。技术宅没有别的本事,算是通过代码来表达一份爱国情怀吧。代码如下:
代码语言:javascript复制# modify ddf and add Taiwan to China
# Taiwan always belong to China
tw = ddf[ddf$Country == 'China',]
tw[1,1] = "Taiwan"
ddf = rbind(ddf,tw)
rownames(ddf) = NULL
# a simple frame of world map
data(wrld_simpl)
plot(wrld_simpl)
读取坐标数据之后呢,这时候数据全部存储在wrld_simpl这个对象中,通过plot()函数可以直接绘图,图形如下,不过稍微简陋了一些。
LEVEL 2
上面的图形看上去稍微有些简陋,而且,还存在一个比较严重的问题。什么问题呢?大家看看本文最上面的NEJM和Lancet杂志中发表的世界地图,和刚才我们绘制出来的世界地图,有什么不一样?提醒一下大家,世界地图中包含几个大洲?
细心地同学可能发现了,在杂志中发表的世界地图是没有南极洲的,为什么呢?第一,南极洲实在是人数稀少,除了企鹅多之外,我们去统计南极洲的某种疾病的发病率,一点意义都没有。第二,南极洲地域广袤,在我们的图形中占据了非常大的版面,使得我们真正关心的地域被压缩。因此,我们在杂志中发表世界地图之前,我们需要做一个事情:把南极洲去除。怎么做呢?
代码语言:javascript复制# this lets us use the contry name vs 3-letter ISO
wrld_simpl@data$id <- wrld_simpl@data$NAME
wrld <- fortify(wrld_simpl, region="id")
wrld <- subset(wrld, id != "Antarctica") # we don't rly need Antarctica
## base map
ggplot()
geom_map(data=wrld, map=wrld, aes(map_id=id, x=long, y=lat), fill="white", color="#7f7f7f", size=0.25)
然后,我们通过ggplot()去做图,放弃base作图系统,这样方便我们做后续的图像调整,结果如下。横纵坐标分别是什么呢?当然是经度和纬度喽。
LEVEL 3
流行病学中绘制地图主要的目的是为了描述发病率一类的信息,然后用不同的颜色将其Highlight出来。我们已经试着把世界地图给绘制出来了,但是还远远没有达到我们的目的:把不同国家标注出来。这个应该怎么操作呢?由于全世界国家众多,为了让大家方便操作,我们罗列一个Excel表格给大家。在表格中有两列数据,第一列包含了240多个国家,按照字母顺序进行排列,第二列就是该国家对应的某个数据值。我们这里的测试数据是随机生成的,没有任何意义,所以等下我们画出来的图可能有点五颜六色,大家做好心理准备。
代码语言:javascript复制# read in data and arguments
ddf = read.xlsx("contries.names.xlsx",sheet = 1)
# 划分数值间隔
max.value = 100 # 最大值
min.value = 0 # 最小值
interval = 20 # 间隔
# add our colored regions
ggplot()
geom_map(data=wrld, map=wrld, aes(map_id=id, x=long, y=lat), fill="white", color="#7f7f7f", size=0.25)
geom_map(data=ddf, map=wrld, aes(map_id=Country, fill=brk), color="white", size=0.25)
有了这个数据之后,我们按照上面的代码,就可以把表格中的数据读入,分割成不同的数值高低,进行绘图了。在这里需要大家在代码中手动做一件事情,由于不知道各位绘图的数据到底长啥样,所以,当你用自己的数据绘图时,请在代码中务必指定一下:最小值,最大值和间隔。我们的测试数据中,最小值0,最大值100,间隔10,当然我们在代码中也都给大家标注了需要大家修改的地方。
说实话,这个图做出来是不是很惊艳?ps,放大起来看,在技术调整过代码后,中国台湾省和咱们中国大陆的颜色完全一致哦。
LEVEL 4
那有没有更加高级的方法,使得数值不要按照分层去进行绘图,而是直接作为一个连续性变量进行绘图呢?这样可能看上去颜色会更流畅一些。当然可以!或者,有没有可能我们自己定义图中的颜色,而不要采用ggplot2已经预定义好的颜色呢?当然也可以喽。上代码:
代码语言:javascript复制# define your own color panel
ggplot()
geom_map(data=wrld, map=wrld, aes(map_id=id, x=long, y=lat), fill="white", color="#7f7f7f", size=0.25)
geom_map(data=ddf, map=wrld, aes(map_id=Country, fill=Value), color="white", size=0.25)
scale_fill_gradientn(colours=colorpanel(75,
low="darkgreen",
mid="yellow",
high="red"))
在代码中,我们完全可以自己定义渐变的颜色,从 "darkgreen"到 “yellow”再到 “red”,这些颜色都是我们可以自己定义的哦。我们这里用的绿到红的颜色渐变,画出来图是这样的:
又一次惊艳了,有木有!这么漂亮的图,画起来不过1-2分钟的时间,简单又实惠哦。
LEVEL 5
那么,经历过这么多次的华丽转变,我们的世界地图已经非常完美了,完全可以用来放到杂志中发表。不过,作为力臻完美的科研猫技术宅,有没有可能更进一步呢?试着把NEJM当中的配色或者Lancet当中的配色直接用到图中去。我们借用了ggsci包中的配色,用到我们的绘图中去,直接上代码:
代码语言:javascript复制# NEJM
ggplot()
geom_map(data=wrld, map=wrld, aes(map_id=id, x=long, y=lat), fill="white", color="#7f7f7f", size=0.25)
geom_map(data=ddf, map=wrld, aes(map_id=Country, fill=brk), color="white", size=0.25)
scale_fill_nejm(alpha = 0.6, name="Value")
theme_void()
# lancet
ggplot()
geom_map(data=wrld, map=wrld, aes(map_id=id, x=long, y=lat), fill="white", color="#7f7f7f", size=0.25)
geom_map(data=ddf, map=wrld, aes(map_id=Country, fill=brk), color="white", size=0.25)
scale_fill_lancet(alpha = 0.6, name="Value")
theme_void()
那么画出来的图是啥样呢?
好了,关于世界地图的绘制,我们先讲到这里。