【问题标题】:Overlaying ggplot, geom_polygon choropleth on ggmap giving errors在ggmap上覆盖ggplot,geom_polygon choropleth给出错误
【发布时间】:2016-10-25 23:44:44
【问题描述】:

我正在努力使用 ggmap、ggplot 和 geom_polygon 在地图上叠加等值线。我对 R 很陌生,所以我不知道如何继续。我正在使用 R Studio。我已经上传了我的数据(这只是假测试数据)、shapefile 和 R 代码here

我有一个阿尔伯塔省人口普查传播区域内的受试者数量数据文件(实际上只是用于测试的随机数)。我也有这些区域的 shapefile。我只对埃德蒙顿市地区感兴趣,如果没有基本的路线图,很难看出我们正在关注城市的哪个部分。

我的代码如下:

#################################
#enter long and lat for appx area of City of Edmonton (this is so we can change it easier later to narrow down if we want)
lat1 <- 53.65 
lat2 <- 53.44 
long1 <- -113.68 
long2 <- -113.35
#################################

xcoords <- c(long2, long1) 
ycoords <- c(lat2, lat1) 

#load data 
data <-read.csv("C:.../fakedata.csv") #data file location 
names(data) 
data

#Load in the shape file 
shp=readShapeSpatial('C:.../gda_048b06a_e.shp') #shape file location

plot(shp) 

这(如预期)给出了该省的所有 DA。但是,这太大了,而且有太多的 DA 让我手动查看哪些在埃德蒙顿,哪些不在,所以我决定根据主题数对它们进行着色,并限制为特别长和纬度:

#Making the data types match 
data$amount <-as.numeric(as.character(data$amount))

is.factor(shp$DAUID) 
data$DAUID<-as.factor(data$DAUID)
is.factor(data$DAUID) #should be True 
is.factor(shp$DAUID) #should be True

shp@data <- left_join(shp@data, data) 
head(shp@data) 
map2 <- fortify(shp, region = "DAUID")

#wait a bit for the fortify to work!
map2 <- rename(map2, DAUID = id) 
map2 <- left_join(map2, shp@data) 
ggplot2 <- ggplot(map2)

#DA map 
datamap <- ggplot2 + geom_polygon(aes(long, lat, group = group, fill = amount),color = "grey") +
scale_fill_gradient(low='white', high='red') + coord_map(xlim =
xcoords,ylim = ycoords) 
datamap

datamap 也是我期望的样子。

#map of Edm 
edmonton <- get_map(location = 'edmonton',maptype="road") 
ggmap(edmonton)

我想做的是将我的 DA 地图叠加到路线图上,这样我们就可以更好地了解暗区和亮区的位置。正是在这一点上,我开始收到错误。

#stuck on how to overlay datamap on top of ggmap(edmonton)

datamap2 <- ggmap(edmonton) + datamap 
datamap2

我收到此错误:

p + o 中的错误:二元运算符的非数字参数另外: 警告消息:不兼容的方法(“+.gg”、“Ops.data.frame”) “+”

我已经用它玩了一个星期,试图将我看到的类似问题的解决方案结合起来,但我似乎无法让它发挥作用。由于我的 shapefile 比我的城市地图大得多,这可能是问题所在吗?我以前从未在 R 中做过这样的事情,而且 GIS 超出了我的技能范围。

我的问题与previous question 类似,但其中的解决方案对我不起作用。如果内容与其他问题过于接近,我深表歉意,并感谢任何指导!

(如果这太长了,我也很抱歉!我不想遗漏任何有用的信息)

【问题讨论】:

    标签: r ggplot2 ggmap


    【解决方案1】:

    我调查了你的案子,最终以我的方式完成了工作。我希望你不要介意。您需要的是在多边形数据的假数据集中有一个带有amount 的数据集。我想这就是您可能使用left_join 的原因。无论如何,每个多边形都有不同数量的数据点。我计算了每个多边形的数据点总数并创建了num。使用它,您可以使用DAUID 创建一个向量,并将id 替换为mymap2。整理好数据后,您将获得带有ggmap 的栅格地图。然后,您在栅格地图上绘制多边形。您可以控制 alpha 值并查看颜色是如何产生的。例如,如果将该值降低到 0.1,则很难看出颜色之间的差异。您可能需要考虑另一种显示道路的方式。一种方法是使用加拿大道路数据集并在埃德蒙顿绘制一些主要道路。希望这会对你有所帮助。

    library(readr)
    library(dplyr)
    library(ggplot2)
    library(ggmap)
    library(rgdal)
    
    #load data
    mydata <- read_csv("fakedata.csv") %>%
              rename(id = DAUID)
    
    #Load in the shape file
    mymap <- readOGR(dsn = ".", layer = "gda_048b06a_e")
    
    mymap@data$DAUID <- as.numeric(as.character(mymap@data$DAUID))
    
    mymap2 <- fortify(mymap)
    
    # Get the number of data points for each polygon
    group_by(mymap2, as.numeric(id)) %>%
    summarize(total = n()) -> num
    
    # Replace id with DAUID, and merge with mydata
    mymap2 %>%
    mutate(id = rep(unique(mymap$DAUID), num$total)) %>%
    left_join(mydata, by = "id") -> mymap2
    
    
    #map of Edm
    edmonton <- get_map(location = "edmonton",maptype = "road")
    
    ggmap(edmonton) + 
    geom_map(data = mymap2, map = mymap2,
             aes(x = long, y = lat, group = group, map_id = id, fill = amount),
             color = "black", size = 0.2, alpha = 0.3) +
    scale_fill_gradient(low = "white", high = "red") +
    coord_map(xlim = c(-113.35, -113.68), ylim = c(53.44, 53.65))
    

    【讨论】:

    • 太完美了!正是我想要做的。我没有意识到问题出在多边形上。非常感谢您的帮助! :)
    • @A.Mullins 很高兴为您提供帮助。我不确定你在做什么。所以我用我自己的方式把它搞定了。我只是认为错误消息表明您尝试添加到 ggmap() 的内容不兼容。如果你看到我的代码,我不会打电话给ggplot(),对吧?在您的情况下,您调用了ggmap(),并尝试添加一个包含ggplot() 的对象。我认为这是错误消息的原因。但最重要的是,您现在拥有了您想要的地图。 :)
    • 供您参考,我正在开发一个包,该包使用单个函数 geom_boundary( )github.com/GL-Li/ggtiger 在 ggmap 上绘制人口普查边界。希望对你有用。
    猜你喜欢
    • 2016-01-14
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2011-08-04
    • 2020-09-23
    • 2015-06-11
    • 1970-01-01
    相关资源
    最近更新 更多