问题描述
我有一个数据框（请参阅下文），具有超过50000个值，每个值与一个位置（纬度，经度）相关联。我想计算5°纬度x 5°经度网格的每个像元的平均值，以创建热图。最终目标是在测深图上绘制此网格。
我看着类似的问题Average values of a point dataset to a grid dataset。但是我无法用自己的数据复制这些示例。遗憾的是，我停留在创建网格的第一步。
我的数据如下：
library(sp)
library(proj4)
coordinates(data) <- c("lon","lat")
proj4string(data) <- CRS("+init=epsg:4326") #defined CRS to WGS 84
df<- data.frame(data)
> head(df)
lon lat value
1 -48.1673562 57.71791 822.9
2 -48.7430053 57.83568 1302.3
3 -48.5662663 57.82087 1508.0
4 -48.3252052 58.29815 224.0
5 -47.1716772 58.42417 38.0
6 -46.4098311 58.67651 431.2
7 -45.8071218 58.70022 365.6
8 -45.5558936 58.46975 50.0
理想情况下，我想使用ggplot2从marmap包中在地图上绘制网格（见下文）：
library(marmap)
library(ggplot2)
atlantic <- getNOAA.bathy(-80,40,90,resolution = 25,keep = TRUE)
atl.df <- fortify(atlantic)
map <- ggplot(atl.df,aes(x=x,y=y)) +
geom_raster(aes(fill=z),data=atl.df) +
geom_contour(aes(z=z),breaks=0,#contour for continent
colour="black",size=1) +
scale_fill_gradientn(values = scales::rescale(c(-5000,1,2400)),colors = c("steelblue4","#C7E0FF","gray40","white"))
解决方法
听起来您想将数值变量（经度和纬度）切成均匀的间隔并汇总每个间隔内的值。以下内容对您有用吗？
library(dplyr)
df2 <- df %>%
mutate(lon.group = cut(lon,breaks = seq(floor(min(df$lon)),ceiling(max(df$lon)),by = 5),labels = seq(floor(min(df$lon)) + 2.5,by = 5)),lat.group = cut(lat,breaks = seq(floor(min(df$lat)),ceiling(max(df$lat)),labels = seq(floor(min(df$lat)) + 2.5,by = 5))) %>%
group_by(lon.group,lat.group) %>%
summarise(value = mean(value),.groups = "drop") %>%
mutate(across(where(is.factor),~as.numeric(as.character(.x))))
样本数据：
set.seed(444)
n <- 10000
df <- data.frame(lon = runif(n,min = -100,max = -50),lat = runif(n,min = 30,max = 80),value = runif(n,min = 0,max = 1000))
> summary(df)
lon lat value
Min. :-99.99 Min. :30.00 Min. : 0.1136
1st Qu.:-87.55 1st Qu.:42.45 1st Qu.: 247.2377
Median :-75.29 Median :55.11 Median : 501.4165
Mean :-75.12 Mean :55.01 Mean : 499.5385
3rd Qu.:-62.69 3rd Qu.:67.63 3rd Qu.: 748.8834
Max. :-50.01 Max. :80.00 Max. : 999.9600
前后数据的比较：
gridExtra::grid.arrange(
ggplot(df,aes(x = lon,y = lat,colour = value)) +
geom_point() +
ggtitle("Original points"),ggplot(df2,aes(x = lon.group,y = lat.group,fill = value)) +
geom_raster() +
ggtitle("Summarised grid"),nrow = 1
)
总是（几乎！），有一个功能。我相信marmap::griddify()是您要寻找的。帮助文件指出：
将不规则间隔的xyz数据转换为适合于创建具有规则间隔的经度和纬度的深水对象的栅格对象。
以下是使用您的坐标的脚本：
library(marmap)
library(ggplot2)
# Create fake data
set.seed(42)
n <- 10000
data_irregular <- data.frame(lon = runif(n,min = -80,max = 40),max = 90),max = 1000))
# Fit data into a grid of 30 cells in longitude and 50 cells in latitude
data_grid <- as.bathy(griddify(data_irregular,nlon = 30,nlat = 50))
fortified_grid <- fortify(data_grid)
# Get bathymetric data to plot continent contours
atlantic <- getNOAA.bathy(-80,40,90,resolution = 25)
atl_df <- fortify(atlantic)
# Plot with ggplot with gridded data as tiles
map <- ggplot(atl_df,aes(x = x,y = y)) +
geom_raster(data = fortified_grid,aes(fill = z)) +
geom_contour(data = atl_df,aes(z = z),breaks = 0,# contour for continent
colour = "black",size = 1) +
scale_fill_gradientn(values = scales::rescale(c(-5000,1,2400)),colors = c("steelblue4","#C7E0FF","gray40","white")) +
labs(x = "Longitude",y = "Latitude",fill = "Value")
map +
theme_bw()
结果如下：
