计算 R 中栅格图层上的点数

Counting number of points on a raster layer in R

我有一张地图,上面有一定数量的点。我想 (1) 计算落在栅格层内的点数,以及 (2) 将这些点提取到数据框。

这是我所做的:

# Packages
library(raster)
library(ggplot2)
library(maptools)
library(tidyverse)
library(dplyr)
library(sp)


# Transform tree ring kml to dataframe
zz<-getKMLcoordinates('treering.kml', ignoreAltitude=TRUE)
l<-as.data.frame(zz)
l<-t(l)

tree <-SpatialPointsDataFrame(l, l,
                                  proj4string = CRS(" +proj=longlat +ellps=WGS84 +datum=WGS84 
                                  +no_defs +towgs84=0,0,0"))


# Get world map
data(wrld_simpl)

# Transform World to raster
r <- raster(wrld_simpl, res = 1)
wrld_r <- rasterize(wrld_simpl, r)

# Import permafrost layer to raster
dist1<-raster("PZI.flt")


# Set CRS
dist1 <- projectRaster(from = dist1, crs = CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs 
+towgs84=0,0,0"))

# Change colours

micolor <- rev(rainbow(12, alpha = 0.35))
transp <- rainbow(12, alpha = 0)
micolor[1:3] <- transp[1]


# Plot all
plot(wrld_r, col = "lightgrey")
plot(dist1, add=TRUE, legend = F, col = micolor)
plot(tree, add=T, pch = 20, col='black', cex=0.2)

我想计算并提取位于这张地图彩色部分的黑点

首先 raster::projectRaster 不是 "set" 投影,而是在给定变换和重采样的情况下重新投影栅格。鉴于此的计算要求,使用 sp::spTransform 重新投影点数据要快得多。一旦您的数据处于同一投影 space,您可以使用 raster::extract 来提取栅格值。栅格外部或无数据 (NA) 区域中的值将被分配 NA 值。您可以使用带有 which.

的简单 NA 索引来删除这些观察结果

看起来您的数据在永久冻土之外可能具有恒定值。一旦确定了这个值是什么(例如 0),您也可以删除这些点。这是一个有效的例子。首先我们添加包并创建一些与您的相似的示例数据。

library(sp)
library(raster)

dist1 <- raster(nrow=20, ncol=20)
  dist1[] <- sample(1:10, ncell(dist1), replace=TRUE)
    dist1[200:400] <- 0

trees <- sampleRandom(dist1, 100, sp=TRUE)

plot(dist1)
  plot(trees,pch=20,col="red",add=TRUE)

现在,我们提取栅格值并查看点对象的维度(请注意,我不必在 raster::extract 函数中使用 sp=TRUE 参数)。

trees@data <- data.frame(trees@data, dist1 = extract(dist1, trees))
  dim(trees)

现在我们创建一个行索引,指示哪些行包含零,确保我们已经识别出行(使用 if 语句),然后删除它们。再次查看对象维度,我们可以看到从原始点数据中删除了多少点。

( idx <- which(trees$dist1 %in% 0) )
  if(length(idx) > 0) trees <- trees[-idx,]
    dim(trees)