4

我正在尝试绘制地图,但我无法弄清楚为什么以下内容不起作用:

这是一个最小的例子

testdf <- structure(list(x = c(48.97, 44.22, 44.99, 48.87, 43.82, 43.16, 38.96, 38.49, 44.98, 43.9), y = c(-119.7, -113.7, -109.3, -120.6,  -109.6, -121.2, -114.2, -118.9, -109.7, -114.1), z = c(0.001216,  0.001631, 0.001801, 0.002081, 0.002158, 0.002265, 0.002298, 0.002334, 0.002349, 0.00249)), .Names = c("x", "y", "z"), row.names = c(NA, 10L), class = "data.frame")

这适用于 1-8 行:

ggplot(data = testdf[1,], aes(x,y,fill = z)) + geom_tile()
ggplot(data = testdf[1:8,], aes(x,y,fill = z)) + geom_tile()

但不适用于 9 行:

ggplot(data = testdf[1:9,], aes(x,y,fill = z)) + geom_tile()

最终,我正在寻找一种在非常规网格上绘制数据的方法。我使用 geom_tile 不是必需的,但是任何点上的空间填充插值都可以。

完整的数据集可作为gist

testdf以上是完整数据集的一小部分,美国的高分辨率栅格(>7500 行)

require(RCurl) # requires libcurl; sudo apt-get install libcurl4-openssl-dev
tmp <- getURL("https://gist.github.com/raw/4635980/f657dcdfab7b951c7b8b921b3a109c7df1697eb8/test.csv")
testdf <- read.csv(textConnection(x))

我试过的:

  1. 使用 geom_point 有效,但没有达到预期的效果:

    ggplot(data = testdf, aes(x,y,color=z)) + geom_point()
    
  2. 如果我将x 或 y 转换为 1:10 的向量,则绘图按预期工作

    newdf <- transform(testdf, y =1:10)
    
    ggplot(data = newdf[1:9,], aes(x,y,fill = z)) + geom_tile()
    
    newdf <- transform(testdf, x =1:10)
    ggplot(data = newdf[1:9,], aes(x,y,fill = z)) + geom_tile()
    

sessionInfo()R version 2.15.2 (2012-10-26) Platform: x86_64-pc-linux-gnu (64-bit)


> attached base packages: [1] stats     graphics  grDevices utils    
> datasets  methods   base     

> other attached packages: [1] reshape2_1.2.2 maps_2.3-0    
> betymaps_1.0   ggmap_2.2      ggplot2_0.9.3 

> loaded via a namespace (and not attached):  [1] colorspace_1.2-0   
> dichromat_1.2-4     digest_0.6.1        grid_2.15.2        
> gtable_0.1.2        labeling_0.1         [7] MASS_7.3-23        
> munsell_0.4         plyr_1.8            png_0.1-4          
> proto_0.3-10        RColorBrewer_1.0-5  [13] RgoogleMaps_1.2.0.2
> rjson_0.2.12        scales_0.2.3        stringr_0.6.2      
> tools_2.15.2
4

4 回答 4

11

您不能使用的原因geom_tile()(或更合适geom_raster()的原因是因为这两个geoms依赖于您的图块是均匀分布的,而事实并非如此。您需要将数据强制转换为点,然后将它们重新采样为均匀分布的栅格,您可以然后用 绘图geom_raster()。您将不得不接受您需要稍微重新采样原始数据才能按照您的意愿绘制它。

您还应该阅读有关地图投影的更多信息raster:::projectionrgdal:::spTransform

require( RCurl )
require( raster )
require( sp )
require( ggplot2 )
tmp <- getURL("https://gist.github.com/geophtwombly/4635980/raw/f657dcdfab7b951c7b8b921b3a109c7df1697eb8/test.csv")
testdf <- read.csv(textConnection(tmp))
spdf <- SpatialPointsDataFrame( data.frame( x = testdf$y , y = testdf$x ) , data = data.frame( z = testdf$z ) )

# Plotting the points reveals the unevenly spaced nature of the points
spplot(spdf)

在此处输入图像描述

# You can see the uneven nature of the data even better here via the moire pattern
plot(spdf)

在此处输入图像描述

# Make an evenly spaced raster, the same extent as original data
e <- extent( spdf )

# Determine ratio between x and y dimensions
ratio <- ( e@xmax - e@xmin ) / ( e@ymax - e@ymin )

# Create template raster to sample to
r <- raster( nrows = 56 , ncols = floor( 56 * ratio ) , ext = extent(spdf) )
rf <- rasterize( spdf , r , field = "z" , fun = mean )

# Attributes of our new raster (# cells quite close to original data)
rf
class       : RasterLayer 
dimensions  : 56, 135, 7560  (nrow, ncol, ncell)
resolution  : 0.424932, 0.4248191  (x, y)
extent      : -124.5008, -67.13498, 25.21298, 49.00285  (xmin, xmax, ymin, ymax)

# We can then plot this using `geom_tile()` or `geom_raster()`
rdf <- data.frame( rasterToPoints( rf ) )    
ggplot( NULL ) + geom_raster( data = rdf , aes( x , y , fill = layer ) )

在此处输入图像描述

# And as the OP asked for geom_tile, this would be...
ggplot( NULL ) + geom_tile( data = rdf , aes( x , y , fill = layer ) , colour = "white" )

在此处输入图像描述

当然我要补充一点,这个数据毫无意义。您真正必须做的是获取 SpatialPointsDataFrame,为其分配正确的投影信息,然后通过 spTransform 转换为 latlong 坐标,然后对转换后的点进行栅格化。确实,您需要了解有关栅格数据的更多信息。你在这里得到的是一个近似值,但最终它并不是数据的真实反映。

于 2013-03-12T00:30:44.050 回答
9

这不是geom_tile()问题的答案,而是另一种绘制数据的方法。

由于您有 30 公里网格的 x 和 y 坐标(我假设在该网格的中间),因此您可以使用geom_point()和绘制数据。您应该选择适当shape=的值。形状 15 将绘制矩形。

另一个问题是 x 和 y 值 - 在绘制数据时,它们应该被绘制为x=yy=x对应于纬度和经度。

coord_equal()将确保有正确的纵横比(我在网上找到了这个以比例为例的解决方案)。

ggplot(data = testdf, aes(y,x,colour=z)) + geom_point(shape=15)+
  coord_equal(ratio=1/cos(mean(testdf$x)*pi/180))

在此处输入图像描述

于 2013-03-11T17:56:53.370 回答
4

回答:

数据被绘制但非常小。


从这里:

"Tile plot as densely as possible, assuming that every tile is the same size.

考虑这个情节

ggplot(data = testdf[1:2,], aes(x,y,fill = z)) + geom_tile()

在此处输入图像描述

上图中有两个图块。geom_tile考虑到每个图块的大小相同,正在尝试使图尽可能密集。在这里,我们可以制作两个这么大的瓷砖而不会重叠。为 4 块瓷砖腾出足够的空间。

试试下面的图,看看结果图告诉你什么:

df1 <- data.frame(x=c(1:3),y=(1:3))
#     df1
#  x   y
#1 1   1
#2 2   2
#3 3   3
ggplot(data = df1[1,], aes(x,y)) + geom_tile()   
ggplot(data = df1[1:2,], aes(x,y)) + geom_tile() 
ggplot(data = df1[1:3,], aes(x,y)) + geom_tile()

与此示例进行比较:

 df2 <- data.frame(x=c(1:3),y=c(1,20,300))
 df2
 # x   y
#1 1   1
#2 2  20
#3 3 300

 ggplot(data = df2[1,], aes(x,y)) + geom_tile()
 ggplot(data = df2[1:2,], aes(x,y)) + geom_tile()
 ggplot(data = df2[1:3,], aes(x,y)) + geom_tile()

请注意,前两个图是相同的df1df2但第三个图df2是不同的。这是因为我们可以制作的最大瓷砖在 (x[1],y[1])和 (之间x[2],y[2])。任何更多它们会重叠,这在这两个瓷砖和最后第三个瓷砖之间留下了很多空间y=300

尽管我不确定这在这里有多合理,但还有一个width参数。geom_tile你确定你不喜欢这种稀疏数据的另一种选择吗?

(您的完整数据仍在绘制中:请参阅ggplot(data = testdf, aes(x,y)) + geom_tile(width=1000)

于 2013-01-24T23:00:37.427 回答
1

如果你想使用 geom_tile 我认为你需要先聚合:

# NOTE: tmp.csv downloaded from https://gist.github.com/geophtwombly/4635980/raw/f657dcdfab7b951c7b8b921b3a109c7df1697eb8/test.csv
testdf <- read.csv("~/Desktop/tmp.csv") 

# combine x,y coordinates by rounding
testdf$x2 <- round(testdf$x, digits=0)
testdf$y2 <- round(testdf$y, digits=0)

# aggregate on combined coordinates
library(plyr)
testdf <- ddply(testdf, c("x2", "y2"), summarize,
                z = mean(z))

# plot aggregated data using geom_tile
ggplot(data = testdf, aes(y2,x2,fill=z)) +
  geom_tile() +
  coord_equal(ratio=1/cos(mean(testdf$x2)*pi/180)) # copied from @Didzis Elferts answer--nice!

一旦我们完成了所有这些,我们可能会得出结论 geom_point() 更好,正如@Didzis Elferts 所建议的那样。

于 2013-03-11T19:21:27.723 回答