Z - R 中多边形（shapefile）的值答案

【问题标题】：Z - Values for polygon (shapefile) in RZ - R 中多边形（shapefile）的值
【发布时间】：2014-01-30 16:12:36
【问题描述】：

我的目标是在 R 中创建一个 3D 可视化。我有一个柏林市区 (Ortsteile) 的 shapefile，并希望将值（居民/平方公里）突出显示为 z 值。我已经将 shapefile 实现到 R 中，并为 desnity（“Einwohnerd”）的值着色，如下所示：

library(rgdal)
library(sp)

berlin=readOGR(dsn="C...etc.", layer="Ortsteile")

berlin@data

col <- rainbow(length(levels(berlin@data$Name)))
spplot(berlin, "Einwohnerd", col.regions=col, main="Ortsteil Berlins", sub="Datensatz der Stadt Berlin", lwd=.8, col="black")

如何将某个多边形（市区）引用到 z 值（居民/km²），如何突出显示这个 z 值？

希望有人能给出答案！最好的问候某人

感谢您的回答，但我仍在努力寻找将密度用作 z 值的最佳方法，以便创建 3D 模型。我发现不可能使用形状的多边形，但可以栅格化多边形并使用矩阵进行不同的透视和旋转。

这是代码，但最终的 3D 可视化看起来不够清晰和不够好。也许最好以另一种方式计算 z 值，这样第一个值就不会开始那么高，或者使用多边形的中心，而不是在 z 方向上绘制一个列：

library(rgdal)
library(sp)

setwd("C:\\...")
berlin=readOGR(dsn="C:\\...\\Ortsteile", layer="Ortsteile") 

col <- rainbow(length(levels(berlin@data$Name)))  
spplot(berlin, "Einwohnerd", col.regions=col, main="Ortsteil Berlins",                 
sub="Datensatz    der Stadt Berlin", lwd=.8, col="black")

library(raster)

raster <- raster(nrows=100, ncols=200, extent(berlin)) 

test <- rasterize(berlin, raster, field="Einwohnerd")

persp(test, theta = 40, phi = 40, col = "gold", border = NA, shade = 0.5)  

for(i in seq(0,90,10)){     
persp(test, theta = 40, phi = i, col = "gold", border = NA, shade = 0.5)
}

library(rgl)         
library(colorRamps)
mat <- matrix(test[], nrow=test@nrows, byrow=TRUE)
image(mat)
persp3d(z = mat, clab = "m")
persp3d(z = mat, col = rainbow(10),border = "black")
persp3d(z = mat, facets = FALSE, curtain = TRUE)

【问题讨论】：

如果您将 shapefile 在线发布到某个地方（例如 Dropbox）并提供链接，您更有可能获得帮助。如果您的人口密度数据位于单独的文件中，您也应该提供该文件。您是在尝试创建 3 维地图（ala Google 地球），还是根据 z 值对多边形进行着色的等值线地图？
这里有一个链接：dropbox.com/sh/r35joqzrq5jnhfp/v_jq90etOM 密度的信息在那里。最后，我尝试使用基于值 Z 着色的多边形创建地图 :-) 谢谢

标签： r 3d visualization

【解决方案1】：

这是你的想法吗？

library(ggplot2)
library(rgdal)           # for readOGR(...) and spTransform(...)
library(RColorBrewer)    # for brewer.pal(...)

setwd("<directory with shapefile>")
map <- readOGR(dsn=".",layer="Ortsteile")
map <- spTransform(map,CRS=CRS("+init=epsg:4839"))
map.data <- data.frame(id=rownames(map@data), map@data)
map.df   <- fortify(map)
map.df   <- merge(map.df,map.data,by="id")
ggplot(map.df, aes(x=long, y=lat, group=group))+
  geom_polygon(aes(fill=Einwohnerd))+
  geom_path(colour="grey")+
  scale_fill_gradientn(colours=rev(brewer.pal(10,"Spectral")))+
  theme(axis.text=element_blank())+
  labs(title="Berlin Ortsteile", x="", y="")+
  coord_fixed()

说明

这是一个很好的问题，因为它提供了一个在 R 中使用 ggplot 的非常基本的等值线图示例。

可以使用readOGR(...) 将Shapefiles 读入R，生成SpatialDataFrame 对象。后者基本上有两个部分：polygons 部分包含多边形边界的坐标，data 部分包含来自 shapefile 中属性表的信息。这些可以分别引用为map@polygons 和map@data。

上面的代码读取 shapefile 并将坐标转换为epsg:4839。然后我们将多边形 ID（存储在行名中）添加到 map@data 中的其他信息中，创建 map.data。然后我们使用ggplot 中的fortify(...) 函数将多边形转换为适合绘图的数据框（map.df）。该数据框有一列 id，对应于 map.data 中的 id 列。然后我们根据id列将属性信息（map.data）合并到map.df中。

ggplot 调用创建地图图层并渲染地图，如下所示：

ggplot:       set the default dataset to map.df; identify x- and y-axis columns
geom_polygon: identify column for fill (color of polygon)
geom_path:    polygon boundaries
theme:        turn off axis text
labs:         title, turn off x- and y-axis labels
coord_fixed:  ensures that the map is not distorted

关于scale_fill_gradientn(...) 的注释：此函数通过插入colours= 参数中提供的调色板来为填充值分配颜色。这里我们使用来自www.colorbrewer.org 的光谱调色板。不幸的是，这个调色板的颜色被尊崇（蓝色 - 红色），所以我们使用 rev(...) 来反转颜色顺序（高=红色，低=蓝色）。如果您更喜欢 matlab 中常见的饱和度更高的颜色，请使用 library(colorRamps) 并将对 scale_fill_gradientn(...) 的调用替换为：

  scale_fill_gradientn(colours=matlab.like(10))+

【讨论】：