R中与ggmap的自定义属性的地理热图

Art*_*sov 6 gis r heatmap ggplot2 ggmap

目标是建立像http://rentheatmap.com/sanfrancisco.html这样的东西

我用ggmap得到了地图,并能够在它上面绘制点.

library('ggmap')
map <- get_map(location=c(lon=20.46667, lat=44.81667), zoom=12, maptype='roadmap', color='bw')
positions <- data.frame(lon=rnorm(100, mean=20.46667, sd=0.05), lat=rnorm(100, mean=44.81667, sd=0.05), price=rnorm(10, mean=1000, sd=300))
ggmap(map) + geom_point(data=positions, mapping=aes(lon, lat)) + stat_density2d(data=positions, mapping=aes(x=lon, y=lat, fill=..level..), geom="polygon", alpha=0.3)
Run Code Online (Sandbox Code Playgroud)

stat_density2d在地图上

这是一个基于密度的漂亮图像.有没有人知道如何制作看起来相同的东西,但是使用position $ property来构建轮廓和比例?

我通过stackoverflow.com彻底查看并没有找到解决方案.

编辑1

positions$price_cuts <- cut(positions$price, breaks=5)
ggmap(map) + stat_density2d(data=positions, mapping=aes(x=lon, y=lat, fill=price_cuts), alpha=0.3, geom="polygon")
Run Code Online (Sandbox Code Playgroud)

五个独立的stat_density图中的结果: 在此输入图像描述

编辑2(来自hrbrmstr)

positions <- data.frame(lon=rnorm(10000, mean=20.46667, sd=0.05), lat=rnorm(10000, mean=44.81667, sd=0.05), price=rnorm(10, mean=1000, sd=300))
positions$price <- ((20.46667 - positions$lon) ^ 2 + (44.81667 - positions$lat) ^ 2) ^ 0.5 * 10000
positions <- data.frame(lon=rnorm(10000, mean=20.46667, sd=0.05), lat=rnorm(10000, mean=44.81667, sd=0.05))
positions$price <- ((20.46667 - positions$lon) ^ 2 + (44.81667 - positions$lat) ^ 2) ^ 0.5 * 10000
positions <- subset(positions, price < 1000)
positions$price_cuts <- cut(positions$price, breaks=5)
ggmap(map) + geom_hex(data=positions, aes(fill=price_cuts), alpha=0.3)
Run Code Online (Sandbox Code Playgroud)

结果是: 在此输入图像描述

它也可以为真实数据创建一个合适的图像.到目前为止,这是最好的结果.欢迎提出更多建议.

编辑3: 这是测试数据和上述方法的结果:

https://raw.githubusercontent.com/artem-fedosov/share/master/kernel_smoothing_ggplot.csv

test<-read.csv('test.csv')
ggplot(data=test, aes(lon, lat, fill=price_cuts)) + stat_bin2d(, alpha=0.7) + geom_point() + scale_fill_brewer(palette="Blues")
Run Code Online (Sandbox Code Playgroud)

在此输入图像描述

我相信应该使用除密度内核之外的其他方法来计算适当的多边形.似乎该功能应该在开箱即用的ggplot中,但我找不到它.

编辑4:我感谢您花时间和精力找出这个看似不太复杂的问题的正确解决方案.我把你的答案投票作为目标的一个很好的近似值.

我发现了一个问题:带圆圈的数据太过于人为,而且这些方法在读取世界数据上表现不佳.

保罗的方法给了我情节: 在此输入图像描述

它似乎捕获了很酷的数据模式.

jazzurro的认可给了我这个情节: 在此输入图像描述

它也得到了模式.但是,这两个图似乎都不如默认的stat_density2d图那么漂亮.我还会等几天看看是否会出现其他解决方案.如果没有,我会将赏金奖给jazzurro,因为这将是我坚持使用的结果.

有一个开放的python + google_maps版本的必需代码.可能有人会在这里找到灵感:https: //github.com/jeffkaufman/apartment_prices

Pau*_*lar 1

在我看来,您所附链接中的地图是使用插值生成的。考虑到这一点,我想知道是否可以通过将插值栅格覆盖到 ggmap 上来实现类似的禁欲主义。

library(ggmap)
library(akima) 
library(raster) 

## data set-up from question
map <- get_map(location=c(lon=20.46667, lat=44.81667), zoom=12, maptype='roadmap', color='bw')
positions <- data.frame(lon=rnorm(10000, mean=20.46667, sd=0.05), lat=rnorm(10000, mean=44.81667, sd=0.05), price=rnorm(10, mean=1000, sd=300))
positions$price <- ((20.46667 - positions$lon) ^ 2 + (44.81667 - positions$lat) ^ 2) ^ 0.5 * 10000
positions <- data.frame(lon=rnorm(10000, mean=20.46667, sd=0.05), lat=rnorm(10000, mean=44.81667, sd=0.05))
positions$price <- ((20.46667 - positions$lon) ^ 2 + (44.81667 - positions$lat) ^ 2) ^ 0.5 * 10000
positions <- subset(positions, price < 1000)

## interpolate values using akima package and convert to raster
r <- interp(positions$lon, positions$lat, positions$price, 
            xo=seq(min(positions$lon), max(positions$lon), length=100),
            yo=seq(min(positions$lat), max(positions$lat), length=100))
r <- cut(raster(r), breaks=5) 

## plot
ggmap(map) + inset_raster(r, extent(r)@xmin, extent(r)@xmax, extent(r)@ymin, extent(r)@ymax) +
  geom_point(data=positions, mapping=aes(lon, lat), alpha=0.2) 
Run Code Online (Sandbox Code Playgroud)

https://i.stack.imgur.com/qzqfu.png

不幸的是,我不知道如何使用 inset_raster 更改颜色或 alpha...可能是因为我对 ggmap 不熟悉。

编辑1

这是一个非常有趣的问题,让我摸不着头脑。当应用于现实世界数据时,插值并没有完全达到我想象的效果;多边形自己接近,jazzurro 看起来肯定好多了!

想知道为什么栅格方法看起来如此锯齿状,我再次查看了您附加的地图,并注意到数据点周围有明显的缓冲区...我想知道是否可以使用一些 rgeos 工具来尝试复制效果:

library(ggmap)
library(raster)
library(rgeos)
library(gplots)

## data set-up from question
dat <- read.csv("clipboard") # load real world data from your link
dat$price_cuts <- NULL
map <- get_map(location=c(lon=median(dat$lon), lat=median(dat$lat)), zoom=12, maptype='roadmap', color='bw')

## use rgeos to add buffer around points
coordinates(dat) <- c("lon","lat")
polys <- gBuffer(dat, byid=TRUE, width=0.005)

## calculate mean price in each circle
polys <- aggregate(dat, polys, FUN=mean)

## rasterize polygons
r <- raster(extent(polys), ncol=200, nrow=200) # define grid
r <- rasterize(polys, r, polys$price, fun=mean) 

## convert raster object to matrix, assign colors and plot
mat <- as.matrix(r)
colmat <- matrix(rich.colors(10, alpha=0.3)[cut(mat, 10)], nrow=nrow(mat), ncol=ncol(mat))
ggmap(map) + 
  inset_raster(colmat, extent(r)@xmin, extent(r)@xmax, extent(r)@ymin, extent(r)@ymax) +
  geom_point(data=data.frame(dat), mapping=aes(lon, lat), alpha=0.1, cex=0.1) 
Run Code Online (Sandbox Code Playgroud)

在此输入图像描述

PS我发现需要将颜色矩阵发送到 inset_raster 来自定义叠加