2011-12-19 112 views
8

我对ggplot比较新,所以请原谅我,如果我的一些问题很简单或根本无法解决。点数据集到网格数据集的平均值

我想要做的是生成一个国家的“热图”,形状的填充是连续的。此外,我的国家形状为.RData。我使用hadley wickham's script将SpatialPolygon数据转换为数据框。我的数据帧的长和经度数据现在看起来像这样

head(my_df) 
long  lat   group 
6.527187 51.87055 0.1 
6.531768 51.87206 0.1 
6.541202 51.87656 0.1 
6.553331 51.88271 0.1 

这个long/lat数据绘制了德国的轮廓。数据框的其余部分在这里省略,因为我认为它不是必需的。对于特定的长/经点,我还有第二个数据框。这看起来像这样

my_fixed_points 
long  lat   value 
12.817  48.917  0.04 
8.533  52.017  0.034 
8.683  50.117  0.02 
7.217  49.483  0.0542 

我想现在要做的,就是颜色根据对所有说谎的该点的一定距离内的固定点的平均值地图上的每一个点。这样我会得到一个(几乎)持续的国家地图的色彩。 我至今是该国的地图绘制与GGPLOT2

ggplot(my_df,aes(long,lat)) + geom_polygon(aes(group=group), fill="white") + 
geom_path(color="white",aes(group=group)) + coord_equal() 

我的第一个想法是产生撒谎已绘制在地图内的点,然后计算该值每产生点my_generated_point像这样

value_vector <- subset(my_fixed_points, 
    spDistsN1(cbind(my_fixed_points$long, my_fixed_points$lat), 
    c(my_generated_point$long, my_generated_point$lat), longlat=TRUE) < 50, 
    select = value) 
point_value <- mean(value_vector) 

我还没有找到一种方法来产生这些点,虽然。就整个问题而言,我甚至不知道是否有可能以这种方式解决问题。我现在的问题是,如果存在一种方法来产生这些点和/或如果有另一种方法来解决。

解决方案

保赐我几乎得到了我想要的。这里有一个荷兰样本数据的例子。

library(ggplot2) 
library(sp) 
library(automap) 
library(rgdal) 
library(scales) 

#get the spatial data for the Netherlands 
con <- url("http://gadm.org/data/rda/NLD_adm0.RData") 
print(load(con)) 
close(con) 

#transform them into the right format for autoKrige 
gadm_t <- spTransform(gadm, CRS=CRS("+proj=merc +ellps=WGS84")) 

#generate some random values that serve as fixed points 
value_points <- spsample(gadm_t, type="stratified", n = 200) 
values <- data.frame(value = rnorm(dim(coordinates(value_points))[1], 0 ,1)) 
value_df <- SpatialPointsDataFrame(value_points, values) 

#generate a grid that can be estimated from the fixed points 
grd = spsample(gadm_t, type = "regular", n = 4000) 
kr <- autoKrige(value~1, value_df, grd) 
dat = as.data.frame(kr$krige_output) 

#draw the generated grid with the underlying map 
ggplot(gadm_t,aes(long,lat)) + geom_polygon(aes(group=group), fill="white") + geom_path(color="white",aes(group=group)) + coord_equal() + 
geom_tile(aes(x = x1, y = x2, fill = var1.pred), data = dat) + scale_fill_continuous(low = "white", high = muted("orange"), name = "value") 

autoKrige Netherlands

+0

请做一个可重现的例子。 – 2011-12-19 15:17:08

+0

我有一种感觉,你正在寻找一种插值算法,请参阅下面我的帖子,以克里格(地质统计学)为例。 – 2011-12-19 15:53:34

+0

很棒,你已经发布了解决方案,+1。我唯一要指出的是它缺少“静音”功能的库(缩放)。 – Eduardo 2013-06-02 19:36:16

回答

15

我想你想要的是沿着这些路线的东西。我预测这种自制软件对于大型数据集的效率会非常低,但它可以用于一个小的示例数据集。我会研究内核密度,也许是raster包。但也许这适合你...

下面的代码片段计算覆盖原始点数据集点的网格镉浓度的平均值。只考虑接近1000米的点。

library(sp) 
library(ggplot2) 
loadMeuse() 

# Generate a grid to sample on 
bb = bbox(meuse) 
grd = spsample(meuse, type = "regular", n = 4000) 
# Come up with mean cadmium value 
# of all points < 1000m. 
mn_value = sapply(1:length(grd), function(pt) { 
    d = spDistsN1(meuse, grd[pt,]) 
    return(mean(meuse[d < 1000,]$cadmium)) 
}) 

# Make a new object 
dat = data.frame(coordinates(grd), mn_value) 
ggplot(aes(x = x1, y = x2, fill = mn_value), data = dat) + 
    geom_tile() + 
    scale_fill_continuous(low = "white", high = muted("blue")) + 
    coord_equal() 

这导致以下图片:

enter image description here

的另一种方法是使用内插算法。克里格就是一个例子。这是使用自动地图包很容易(临场自我宣传:),我写的包):

library(automap) 
kr = autoKrige(cadmium~1, meuse, meuse.grid) 
dat = as.data.frame(kr$krige_output) 

ggplot(aes(x = x, y = y, fill = var1.pred), data = dat) + 
    geom_tile() + 
    scale_fill_continuous(low = "white", high = muted("blue")) + 
    coord_equal() 

这导致了下面的图片:

enter image description here

然而,如果没有知识到这张地图你的目标是什么,我很难看到你想要的东西。

2

slideshow提供了另一种方法 - 有关该方法的描述的page 18和有关幻灯片制造商的结果看起来如何的page 21

但请注意,幻灯片制作工具使用sp包和spplot函数,而不是ggplot2及其绘图函数。

+0

...此外,spplot使用引擎盖下的格子。 – 2011-12-20 10:35:36