Извлечь растр из списка объектов SpatialPolygonsDataFrame в R

Question

Я пытаюсь извлечь суммарные значения растровых ячеек из одного большого файла для различных объектов SpatialPolygonsDataFrames (SPDF) в R, хранящихся в списке, а затем добавить извлеченные значения в таблицы атрибутов объектов SPDF. Я хотел бы повторить этот процесс, и понятия не имею, как это сделать. Я нашел эффективное решение для нескольких полигонов, хранящихся в одном объекте SPDF (см .: https://gis.stackexchange.com/questions/130522/increasing-speed-of-crop-mask-extract-raster-by-many-polygons-in-r),, но не знаю, как применить процедуру crop> mask> extract к СПИСКУ объектов SPDF каждый из которых содержит несколько полигонов. Вот воспроизводимый пример:

library(maptools)  ## For wrld_simpl
library(raster)

## Example SpatialPolygonsDataFrame
data(wrld_simpl) #polygon of world countries
bound <- wrld_simpl[1:25,] #country subset 1  
bound2 <- wrld_simpl[26:36,] #subset 2

## Example RasterLayer
c <- raster(nrow=2e3, ncol=2e3, crs=proj4string(wrld_simpl), xmn=-180, 
xmx=180, ymn=-90, ymx=90)
c[] <- 1:length(c)

#plot, so you can see it
plot(c)    
plot(bound, add=TRUE) 
plot(bound2, add=TRUE, col=3) 

#make list of two SPDF objects
boundl<-list()
boundl[[1]]<-bound1
boundl[[2]]<-bound2

#confirm creation of SPDF list
boundl

Далее я хотел бы запустить весь список в формате forloop. Для одного SPDF из списка, кажется, работают следующие серии функций:

clip1 <- crop(c, extent(boundl[[1]])) #crops the raster to the extent of the polygon, I do this first because it speeds the mask up
clip2 <- mask(clip1, boundl[[1]]) #crops the raster to the polygon boundary 
extract_clip <- extract(clip2, boundl[[1]], fun=sum)  
#add column + extracted raster values to polygon dataframe
boundl[[1]]@data["newcolumn"] = extract_clip

Но когда я пытаюсь изолировать первую функцию для списка SPDF (raster::crop), она не возвращает растровый объект:

crop1 <- crop(c, extent(boundl[[1]])) #correctly returns object class 'RasterLayer'
cropl <- lapply(boundl, crop, c, extent(boundl)) #incorrectly returns objects of class 'SpatialPolygonsDataFrame'

Когда я пытаюсь изолировать функцию маски для списка SPDF (raster::mask), она возвращает ошибку:

maskl <- lapply(boundl, mask, c) 
#Error in (function (classes, fdef, mtable)  : unable to find an inherited method for function ‘mask’ for signature ‘"SpatialPolygonsDataFrame", "RasterLayer"’

Я хотел бы исправить эти ошибки и эффективно выполнить всю процедуру внутри одного цикла (т. Е. crop> mask> extract> добавить извлеченные значения в таблицы атрибутов SPDF. Я действительно новичок в R и не знаю, куда идти отсюда. Пожалуйста, помогите!

Answer 1

Один из подходов состоит в том, чтобы взять то, что работает, и просто поместить нужную строку «обрезать -> маска -> извлечь -> добавить» в цикл for:

for(i in seq_along(boundl)) {
    clip1 <- crop(c, extent(boundl[[i]])) 
    clip2 <- mask(clip1, boundl[[i]])
    extract_clip <- extract(clip2, boundl[[i]], fun=sum)  
    boundl[[i]]@data["newcolumn"] <- extract_clip
}

Можно ускорить цикл при параллельном выполнении, например, с помощью пакета R foreach . И наоборот, прирост скорости при использовании lapply() вместо цикла for будет небольшим.

Почему возникает ошибка:

cropl <- lapply(boundl, crop, c, extent(boundl)) 

применяет функцию crop() к каждому элементу списка boundl. Выполнена операция

tmp <- crop(boundl[[1]], c)
## test if equal to first element
all.equal(cropl[[1]], tmp) 
[1] TRUE

Для получения желаемого результата используйте

cropl <- lapply(boundl, function(x, c) crop(c, extent(x)), c=c)
## test if the first element is as expected
all.equal(cropl[[1]], crop(c, extent(boundl[[1]]))) 
[1] TRUE

Примечание:

Использование c для обозначения R-объекта - предпочтительный выбор, поскольку его легко спутать с c().