Количество точек выборки в ячейке сетки

Question

Рассчитать общее количество точек выборки в каждой ячейке сетки пространственной сетки.

Я хотел бы создать сетку и рассчитать общее количество точек выборки в каждой ячейке сетки. Я создал случайно сгенерированные данные и сетку и попытался вычислить количество записей в ячейках сетки, используя как sf, так и растровые пакеты, используя предыдущие подобные вопросы SO, но без успеха. Я также изучил функцию извлечения. Я довольно новичок в пространственном анализе.

 library(sf)
library(raster)
library(tidyverse)
library(mapview)
library(mapedit)

#Trial with sf package 
# load some spatial data. Administrative Boundary
#https://stackoverflow.com/questions/41787313/how-to-create-a-grid-of-       spatial-points
aut <- getData('GADM', country = 'aut', level = 0)
aut <- st_as_sf(aut)
#Try with polygons
grid <- aut %>% 
 st_make_grid(cellsize = 0.5, what = "polygons") %>% 
  st_intersection(aut)                               

#fake data
 lat<-runif(1000, 46.5, 48.5)
 lon<-runif(1000, 13,16)
pos<-data.frame(lat,lon)

 ggplot() + 
  geom_sf(data = aut) + 
  geom_sf(data = grid)+
geom_point(data=pos, aes(lon, lat)) 
#how to count number of records within each cell?  
 ########################################
#Trial with raster package
#https://stackoverflow.com/questions/32889531/r-how-can-i-count-how-   many-points-are-in-each-cell-of-my-grid
 r<-raster(xmn=13, ymn=46.5, xmx=16, ymx=48.5, res=0.5)
r[] <- 0
#How do I use the pos data here
xy <- spsample(as(extent(r), 'SpatialPolygons'), 100, 'random')
tab <- table(cellFromXY(r, xy))
r[as.numeric(names(tab))] <- tab
plot(r)
points(xy, pch=20)
d <- data.frame(coordinates(r), count=r[])

Я хотел бы получить таблицу с количеством точек отбора проб.

Answer 1

подсчет lengths из st_intersects (примечание: не st_intersection) даст вам вектор точек, содержащихся в каждой ячейке сетки:

library(sf)
library(raster)
library(tidyverse)
library(mapview)
library(mapedit)

#Trial with sf package 
# load some spatial data. Administrative Boundary
#https://stackoverflow.com/questions/41787313/how-to-create-a-grid-of-       spatial-points
aut <- getData('GADM', country = 'aut', level = 0)
aut <- st_as_sf(aut)
#Try with polygons
grid <- aut %>% 
  st_make_grid(cellsize = 0.5, what = "polygons") %>% 
  st_intersection(aut)                               

#fake data
lat<-runif(1000, 46.5, 48.5)
lon<-runif(1000, 13,16)
pos<-data.frame(lat,lon)

pos = st_as_sf(pos, coords = c("lon", "lat"), crs = 4326)

tab = st_intersects(grid, pos)
lengths(tab)
[1]  0  0  0  0  4 24 23 34 23 13 14  0  0  0  0  0  0  0  3 38 40 48 46 47 33  0  0  0  0  0  0  0
[33]  0 35 48 51 35 38 44  0  0  0  0 44 43 41 53 44 32  0  0  0  0  8  8 10 12  7  0  0  0  0  0

Если вы затем захотите связать это ссетка как sf объект, который вы можете сделать:

grid = st_sf(n = lengths(tab), geometry = st_cast(grid, "MULTIPOLYGON"))

mapview(grid, zcol = "n")

Answer 2

Пример данных

library(raster)   
aut <- getData('GADM', country = 'aut', level = 0)
r <- raster(aut, res=0.5)
lat <- runif(1000, 46.5, 48.5)
lon <- runif(1000, 13,16)
# note that you should use (lon, lat), in that order!
pos <- data.frame(lon, lat)

Решение

r <- rasterize(pos, r, fun="count")
plot(r)

Чтобы получить столик, вы можете сделать

x <- rasterToPoints(r)
z <- cbind(cell=cellFromXY(r, x[,1:2]), value=x[,3])
head(z)
#     cell value
#[1,]   22    4
#[2,]   23   45
#[3,]   24   36
#[4,]   25   52
#[5,]   26   35
#[6,]   27   38

Или, альтернативно, na.omit(cbind(1:ncell(r), values(r)))

Answer 3

Попробуйте

ggplot(pos, aes(x = lon, y=lat)) + 
geom_bin2d(binwidth = 2) +
stat_bin_2d(aes(label=stat(count)), binwidth=2, geom="text", position="identity") +
scale_fill_gradient(low = "white", high = "red")