R矢量地图栅格化（将shapefile转换成raster）

Python017

R矢量地图栅格化（将shapefile转换成raster）,第1张

在处理地图数据时候，经常会碰到 shp 与 raster 两种格式。通常r中应用较多的为raster栅格数据。shp文件太大，读取也不方便。逐渐被 GeoJSON 替代，用sf去处理与读取。

R在读取shp时候，处理，或者画图都会碰到，反应迟钝问题。所以，我们有时候会根据需要，将shp文件转成raster，不仅可视化快，还可方便数据处理与提取。shp文件转成raster主要解决以下问题：

下面就来介绍，如何根据shp文件，转成raster及在转换过程中碰到的一些问题。

利用 raster 包自带的数据进行演示。读取的是 SpatialPolygonsDataFrame ，关于如何读取shp文件，可以用rgdal与sf的命令。

关键是 rasterize， rasterize(shape, r, 1) 里面有三个主要参数：

那如果我们需要根据shp里面的地区数来生成不同的value呢，意思就是，不用地区value不一样，不应该是统一值。

有时候生成的raster里面有NA数据，那么如何替换掉呢，(reclassify)[ http://search.r-project.org/library/raster/html/reclassify.html] 可以实现该过程。主要参数cbind（0,a,b）意思是将0-a的数值全部变成b。

具体参见： ?reclassify

下面我么将NA替换成0，或者，value=12的替换成100.

转换成raster最终目的是实现数据的提取。譬如现在有两个点，如何提取对应点上的value。

如果是shp文件，操作比较麻烦点，又是还会提取出NA。转换Raster以后，就更方便了。

上面的图太模糊了，那我们设置res就好。

res精度提高，运行速度会下降，尤其是遇到很大的shp数据时候。

一般R里面加载shp超过50M，系统就会迟钝。

rasterize 里面还可以设置field=1.可以达到同样效果。

假设有两组栅格数据，一组代表2019年中国每月降雨量，一组代表2019年中国每月植被叶面积指数（LAI）。想要得到中国月降水量与LAI的相关性分布，那么需要对两组栅格数据对应的栅格点进行逐栅格的相关性分析。

将降水数据导入栅格栈中，这个过程可以理解为将降水数据按时间顺序从上到下堆叠。同理，按相同的时间顺序将LAI数据堆叠。值得一提的是，stack（）函数在堆叠栅格数据时是按文件名拼音和数字大小顺序自动堆叠的，具体规则可以亲自尝试。最后，将这两个栅格栈合并成一个。

对相关性分析函数稍作改变。

以上方法是可以推广的，线性回归函数lm（）和相关性分析函数cor（）的输入都可以是向量，因此只要函数支持向量输入，理论上讲都可以类比上述过程实现。但是如果函数只支持数据框输入，如gbm包中的函数gbm（），那就只能另辟蹊径了。

不知道你问的是在同一图形中添加点（类似画散点图）还是要把屏幕一分为多。

添加的话，用完plot，添加点用points，添加线用lines。

简单点的应用类似：

plot（X,Y）

points（X，Y1）

这样的感觉。

一分为多的话，用split.screen。

上我自己最近写的代码做个例子：

jpeg(filename="geeseP3.jpeg") #画jpeg图

split.screen(c(1,2)) #分屏幕为左右两边

screen(1) #屏幕1预备输出

plot(X2,Y,type="p",xlab="X2",ylab="Y",main="Plotting of X2 and Y")

screen(2) #屏幕2预备输出

plot(X2,Y,type="l")

dev.off()

画出来的图大概是这个感觉：

不喜欢这个比例的话，也有命令可以调节图片长宽比例。这个略去不提。

另外推荐这个

http://wenku.baidu.com/link?url=mMnroYY14th1qiKzsFnTUVceptBVugQsrLbYFItaqMN25xftBQlMBThtyW5fsmIgkMWcWbkXyozKR85SFEb7VwDUhekqSBVDuOvskifRo7W

里面有列举了一些画图用的函数。