地理 heat/contour 地图空间插值的最佳方法？

Question

我想使用 ggplot2 和 ggmap 之类的东西来 生成任意值的热图 例如属性价格每平方米在街道级别的地理区域（具有高分辨率）。

不幸的是，这项任务似乎相当困难，因为虽然 ggplot2 可以生成很好的密度图，但如果没有事先插值，它似乎无法像这样可视化空间数据。

为此，我使用了库 akima（不规则数据的网格双变量插值）和 mgcv（具有集成平滑度估计的广义加性模型），但是我对插值方法的了解很一般充其量，我能够产生的结果还不够令人满意。

考虑以下示例：

数据

library(ggplot2)
library(ggmap)

## data simulation
set.seed(1945)

df <- tibble(x = rnorm(500, -0.7406, 0.03),
             y = rnorm(500, 51.9976, 0.03),
             z = abs(rnorm(500, 2000, 1000)))

地图、散点图、密度图

## ggmap
map <- get_map("Bletchley Park, Bletchley, Milton Keynes", zoom = 13, source = "stamen", maptype = "toner-background")
q <- ggmap(map, extent = "device", darken = .5)

## scatterplot over map
q + geom_point(aes(x, y), data = df, colour = z)

## classic density heat map
q + 
  stat_density2d(aes(x=x, y=y, fill=..level..), data=df, geom="polygon", alpha = .2) + 
  geom_density_2d(aes(x=x, y=y), data=df, colour = "white", alpha = .4) +
  scale_fill_distiller(palette = "Spectral")

如您所见，所选区域的数据相当密集，密度热图看起来很棒，具有圆边和闭合曲线（除了一些最外层）。

使用 akima 进行插值和绘图

## akima interpolation
library(akima)

df_akima <-interp2xyz(interp(x=df$x, y=df$y, z=df$z, duplicate="mean", linear = T,
                             xo=seq(min(df$x), max(df$x), length=200),
                             yo=seq(min(df$y), max(df$y), length=200)), data.frame=TRUE)

## akima plot
q +
  geom_tile(aes(x = x, y = y, fill = z), data = df_akima, alpha = .4) +
  stat_contour(aes(x = x, y = y, z = z, fill = ..level..), data = df_akima, geom = 'polygon', alpha = .4) +
  geom_contour(aes(x = x, y = y, z = z), data = df_akima, colour = 'white', alpha = .4) +
  scale_fill_distiller(palette = "Spectral", na.value = NA)

这会产生一个密集的插值网格（以确保足够的分辨率），虽然下面的平铺图是可以接受的，但等高线图太参差不齐，许多曲线没有闭合。

使用 linear = F 的非线性插值更平滑，但显然牺牲了分辨率并且随着数字（z 的负值）变得疯狂。

使用 mgcv 进行插值和绘图

## mgcv interpolation
library(mgcv)

gam <- gam(z ~ s(x, y, bs = 'sos'), data = df)
df_mgcv <- data.frame(expand.grid(x = seq(min(df$x), max(df$x), length=200),
                                  y = seq(min(df$y), max(df$y), length=200)))
resp <- predict(gam, df_mgcv, type = "response")
df_mgcv$z <- resp

## mgcv plot
q +
  geom_tile(aes(x = x, y = y, fill = z), data = df_mgcv, alpha = .4) +
  stat_contour(aes(x = x, y = y, z = z, fill = ..level..), data = df_mgcv, geom = 'polygon', alpha = .4) +
  geom_contour(aes(x = x, y = y, z = z), data = df_mgcv, colour = 'white', alpha = .4) +
  scale_fill_distiller(palette = "Spectral", na.value = NA)

使用 mgcv 的相同过程会产生漂亮且平滑的绘图，但分辨率要低得多，而且实际上所有曲线都没有闭合。

问题

能否请您提出更好的方法或修改我的尝试以获得类似于第一个的图（干净、连接、平滑的高分辨率线条）？
是否可以关闭曲线，例如在最后一个图中（阴影区域应该在图像边界之外计算）？

感谢您的宝贵时间！

Answer 1

抱歉，我目前无法运行您的示例提供详细信息。但是请尝试使用 automap 包中的 autoKrige()。

克里金法是一种很好的插值方法。只要确保您的数据符合要求即可。这是一个很好的指南： https://gisgeography.com/kriging-interpolation-prediction/

Answer 2

您的地图的问题不是您使用的插值方法，而是 ggplot 显示密度线的方式。这是对此的回答：.

密度线超出地图范围，因此任何超出绘图区域的多边形都会被不恰当地渲染（ggplot 将使用相应级别的下一个点关闭多边形）。这在您的第一张地图上显示不多，因为插值分辨率很低。

Andrew提出的技巧是先扩大绘图区域，使密度线正确渲染，然后切断显示区域以隐藏多余的space。由于我用你的第一个例子测试了他的解决方案，这里是代码：

q + 
  stat_density2d(
    aes(x = x, y = y, fill = ..level..),
    data = df,
    geom = "polygon",
    alpha = .2,
    color = "white",
    bins = 20
  ) + 
  scale_fill_distiller(
    palette = "Spectral"
  ) +
  xlim(
    min(df$x) - 10^-5,
    max(df$x) + 10^-5
  ) +
  ylim(
    min(df$y) - 10^-3,
    max(df$y) + 10^-3
  ) +
  coord_equal(
    expand = FALSE,
    xlim = c(-.778, -.688),
    ylim = c(51.965, 52.03)
  )

唯一的区别是我使用 min()- / max() + 而不是固定数字和 coord_equal 来确保地图没有失真。此外，我手动指定了更多的级别（使用 bin），因为通过增加绘图区域，stat_density 会自动选择较低的分辨率。

至于最佳插值方法，这取决于您的 objective 和您拥有的数据类型。问题不是地图的最佳方法是什么，而是数据的最佳方法是什么。这是一个非常广泛的问题，超出了本 space 的范围。但这里有一个很好的指南：http://www.rspatial.org/analysis/rst/4-interpolation.html

关于如何使用 ggplot 在 R 中制作好的地图的一般想法：http://spatial.ly/r/

地理 heat/contour 地图空间插值的最佳方法？

Best method of spatial interpolation for geographic heat/contour maps?

interpolation

r

geospatial

ggplot2

ggmap

数据

地图、散点图、密度图

使用 akima 进行插值和绘图

使用 mgcv 进行插值和绘图

问题