使用marmap包和getNOAA.bathy在ggplot中添加水深线

Question

我想要在地图上添加水深线。我正在绘制海岸线附近的点，并且我们对它们距离大陆架有多近很感兴趣。我看过一个叫做 Marmap 的包，但现在我正在使用 ggplot，因为它能提供更高的分辨率。

获取水深线的代码如下：

library(marmap)

Peru.bath <- getNOAA.bathy (lon1 = -90, lon2 = -70, lat1 = -20, 
                            lat2 = -2, resolution = 10) 

plot(Peru.bath)

我正在使用的代码如下，我想要添加海底地形线：

coast_map <- fortify(map("worldHires", fill=TRUE, plot=FALSE))
gg <- ggplot()
gg <- gg + geom_map(data=coast_map, map=coast_map,
                    aes(x=long, y=lat, map_id=region),
                    fill="white", color="black") +
  theme(panel.background = element_blank()) + 
  theme(panel.grid.major = element_blank()) + 
  theme(panel.grid.minor = element_blank()) +
  theme(axis.text.x = element_blank(), axis.text.y = element_blank(), axis.ticks = element_blank())
gg <- gg + xlab("") + ylab("")
gg <- gg + geom_map(data=data.frame(region="Peru"), map=coast_map,
                    aes(map_id=region), fill="gray") 
gg <- gg + xlim(-90,-70) + ylim(-20,-2)
gg <- gg + coord_map()
gg

因此，我认为它会是

gg <- gg + Peru.bath

然而，我遇到了“错误：不知道如何将Peru.bath添加到图形中”的问题。

请注意，我没有水深数据，我只是希望在我创建的地图上绘制已知的海洋架线，如果可能的话。

- CDav

Answer 1

我刚刚在Github上更新了marmap的开发版本。您可以通过以下方式进行安装：

library(devtools)
install_github("ericpante/marmap")

现在包括了使用ggplot2绘制bathy对象的autoplot.bathy()函数。请务必查看其帮助文件和示例，以了解可能性。这是一个使用数据集dat的示例：

library(marmap) ; library(ggplot2)
dat <- getNOAA.bathy(-90,-70,-20,-2,res=4, keep=TRUE)

# Plot bathy object using custom ggplot2 functions
autoplot(dat, geom=c("r", "c"), colour="white", size=0.1) + scale_fill_etopo()

- Benoit

Answer 2

如果您想要（i）在地图上添加等深线和（ii）使用get.depth()确定数据点的深度，最好使用标准绘图方法（marmap已经与标准绘图相兼容）。事实上，“更好的分辨率”与基本绘图或ggplot2无关。在您的例子中，您绘制的海岸线来源于mapdata包中的“worldHires”数据集，与ggplot2无关。确实，在marmap绘图中也可以添加相同的海岸线。

以下是一些代码，可以使用基本绘图和marmap产生两个相当不错的地图：

library(marmap) ; library(mapdata)

# Get bathymetric data
dat <- getNOAA.bathy(-90,-70,-20,-2,res=4, keep=TRUE)

# Create nice color palettes
blues <- c("lightsteelblue4", "lightsteelblue3", "lightsteelblue2", "lightsteelblue1")
greys <- c(grey(0.6), grey(0.93), grey(0.99))

## First option for plotting
plot(dat, land=TRUE, n=100, lwd=0.03)
map("worldHires", res=0, add=TRUE)

# Second option
plot(dat, im=TRUE, land=TRUE, bpal=list(c(min(dat),0,blues),c(0,max(dat),greys)), lwd=.05, las=1 )
map("worldHires", res=0, lwd=0.7, add=TRUE)

# Add -200m and -1000m isobath
plot(dat, deep=-200, shallow=-200, step=0, lwd=0.5, drawlabel=TRUE, add=TRUE)
plot(dat, deep=-1000, shallow=-1000, step=0, lwd=0.3, drawlabel=TRUE, add=TRUE)

请注意，这里使用的分辨率不是最高的。在getNOAA.bathy()函数中，res参数设置为4。这将下载一个2.7Mb的数据集，并可以通过将keep参数设置为TRUE来保存到本地。最高可能的分辨率是res=1，但在我看来，对于如此广泛的地理范围，这已经超过需求。以下代码生成了下面两个图表：

- Benoit

---

这太棒了@Benoit！非常感谢。然而，自从发帖以来，我实际上已经继续使用Google地图作为基础包进行ggplot的努力，因为卫星图像显示了一些接近海底地形的数据-有没有办法将其添加为ggplot基础地图？我尝试过fortify，但它说ggplot不认识bathy类？非常感谢。 - CDav

请看下面的第二个答案：我已经包含了一个使用新函数在ggplot2中绘制bathy对象的示例。 - Benoit

Answer 3

早上时间有点不够，但这里的内容可以帮助你入门（当然，其他R地理专家也可以进一步完善）：

library(maps)
library(mapdata)
library(ggplot2)
library(marmap)
library(Grid2Polygons)

coast_map <- fortify(map("worldHires", fill = TRUE, plot = FALSE))

Peru.bath <- getNOAA.bathy (lon1 = -90, lon2 = -70, lat1 = -20,
                            lat2 = -2, resolution = 10)

peru_bathy_df <- Grid2Polygons(as.SpatialGridDataFrame(Peru.bath), 
                               level=TRUE, pretty=TRUE)
peru_bathy_map <- fortify(peru_bathy_df)

gg <- ggplot()
gg <- gg + geom_map(data=peru_bathy_map, map=peru_bathy_map,
                    aes(map_id=id), color="black", fill="white")
gg <- gg + geom_map(data=coast_map, map=coast_map,
                    aes(x=long, y=lat, map_id=region),
                    fill="white", color="black")
gg <- gg + geom_map(data=data.frame(region="Peru"), map=coast_map,
                    aes(map_id=region), fill="steelblue")
gg <- gg + xlim(-90,-70) + ylim(-20,-2)
gg <- gg + coord_map()
gg <- gg + theme_bw()
gg

显然，你希望得到比这更好的图片，但基本思路是将其转换为ggplot可以处理的对象（即SpatialPolygonsDataFrame）。在bathy和Grid2Polygons帮助文件中有很多非ggplot的好例子。

注意：这些转换/渲染需要一些时间，并且bathy示例展示了一种无需ggplot的方法，速度会快得多。

- hrbrmstr

---

非常感谢这个 - 看起来很棒！但是，我想知道是否有可能保留线条并添加深度信息。我已经使用了get.depth函数来获取我想绘制的每个数据点的深度信息，但不确定如何将其添加到图表中（最好使用ggplot）...也许@Benoit能够指导我正确的方向？ - CDav