使用 extract() 函数在点处提取栅格数据:如何防止 NA 的输出

Extracting Raster Data at Points using the extract() function: How To Prevent the Output of NA's

我有一个 data.frame 具有三个变量(ID、经度和纬度)的海豚观察。

我还有 59 个 AQUA/MODIS 包含海面温度 (SST) 的 netCDF 文件。我想提取海豚观测位置的 SST 栅格数据。

海豚数据

 d <- read.csv("dolphins.csv")
 str(d) 
 #  'data.frame':   650 obs. of  4 variables:
 #$ ID       : num  1 2 3 4 5 6 7 8 9 10 ...
 #$ Date     : Factor w/ 73 levels "1/24/17","1/24/18",..: 67 67 67 8 8 8 8 8 8 8 ...
 #$ Latitude : num  42 42 42 42 42 ...
 #$ Longitude: num  19.1 19.1 19.1 19.1 19.1 ...

SST 栅格

library(terra)
filenames = list.files('Ocean_ColorSST_2016_2021',pattern='*.nc',full.names=TRUE)
SSTs <- rast(filenames, "sst")
SSTs
#class       : SpatRaster 
#dimensions  : 4320, 8640, 59  (nrow, ncol, nlyr)
#resolution  : 0.04166667, 0.04166667  (x, y)
#extent      : -180, 180, -90.00001, 90  (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84 
#sources     : AQUA_MODIS.20160901_20160930.L3m.MO.SST.sst.4km.nc:sst  
#              AQUA_MODIS.20161001_20161031.L3m.MO.SST.sst.4km.nc:sst  
#              AQUA_MODIS.20161101_20161130.L3m.MO.SST.sst.4km.nc:sst  
#              ... and 56 more source(s)
#varnames    : sst (Sea Surface Temperature) 
#              sst (Sea Surface Temperature) 
#              sst (Sea Surface Temperature) 
#              ...
#names       :      sst,      sst,      sst,      sst,      sst,      sst, ... 
#unit        : degree_C, degree_C, degree_C, degree_C, degree_C, degree_C, ... 

创建一个 SpatialPointsDataFrame

library(sp)
points_spdf_M <- d
coordinates(points_spdf_M) <- ~ Latitude + Longitude
crs(points_spdf_M) <- crs(SSTs)
points_spdf_M
#class       : SpatialPoints 
#features    : 650 
#extent      : 41.5978, 42.67778, 18.24337, 19.99933  (xmin, xmax, ymin, ymax)
#crs         : +proj=longlat +datum=WGS84 +no_defs 

从海豚 ID returns NAs

中提取光栅数据
library(raster)
ncin_SST <- stack(SSTs)
Extract_SST_M <- extract(ncin_SST, points_spdf_M)
head(Extract_SST_M)
#sst.1 sst.2 sst.3 sst.4 sst.5 sst.6 sst.7 sst.8 sst.9 sst.10 sst.11 sst.12 sst.13 sst.14 sst.15 sst.16 sst.17 sst.18
#  [1,]    NA    NA    NA    NA    NA    NA    NA    NA    NA     NA     NA     NA     NA     NA     NA     NA     NA     NA
# [2,]    NA    NA    NA    NA    NA    NA    NA    NA    NA     NA     NA     NA     NA     NA     NA     NA     NA     NA
# [3,]    NA    NA    NA    NA    NA    NA    NA    NA    NA     NA     NA     NA     NA     NA     NA     NA     NA     NA

我试过的另一种方法

v <- vect(d, geom=c("Longitude", "Latitude"))
e <- extract(SSTs, v)
head(e)

输出 ?

你像这样创建一个 SpatialPointsDataFrame

coordinates(Final_M_Points) <- ~ Latitude + Longitude

它应该在哪里

coordinates(Final_M_Points) <- ~ Longitude + Latitude

您的工作流程得到简化 --- 仅使用 terra

library(terra)
filenames <- list.files("---")
SSTs <- rast(filenames, "sst")

d <- read.csv("file.csv")
v <- vect(d, geom=c("Longitude", "Latitude"))

e <- extract(SSTs, v)

最好用terra,避免spraster。在这种情况下,您会收到经度和纬度颠倒的警告。