Как найти среднее растровое значение области, определенной шейп-файлом, используя R?

19

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

Как извлечь значения из растров и импортировать их в R так, чтобы я мог их использовать?

nickves
источник

Ответы:

23

Вот пример кода. Довольно просто адаптировать этот код для работы в цикле для обработки всех ваших растров. Если ваши растры имеют общий экстент и разрешение, вы можете создать растровый стек и циклически проходить по полосам в стеке. Чтобы создать вектор, содержащий все растры в каталоге, в определенном формате, вы можете использовать «list.files», а затем передать этот вектор в стек.

Пример:

rlist=list.files(getwd(), pattern="img$", full.names=TRUE) 

r <- stack(rlist)   


    # Add required libraries
    require(raster)
    require(sp)
    require(rgdal)

    # Set working directory, raster, in and out shapefiles
    setwd("C:/test")
    inshp="MyPolys"
    outshp="PolyMeans"
    rdata <- "Year2012.img"

    # Read polygon feature class shapefile
    sdata <- readOGR(dsn=getwd(), layer=inshp)

    # Read raster
    r <- raster(rdata)

    # Extract raster values to list object
    r.vals <- extract(r, sdata)

    # Use list apply to calculate mean for each polygon
    r.mean <- lapply(r.vals, FUN=mean)

    # Join mean values to polygon data
    sdata@data <- data.frame(sdata@data, m2012=r.mean)

    # Write results
    writeOGR(sdata, getwd(), outshp, driver="ESRI Shapefile", check_exists=TRUE, 
         overwrite_layer=TRUE)
Джеффри Эванс
источник
1
+1 - отличное решение, хорошо продуман и полный ответ.
Симбамангу
Если у меня есть мультиполигональный файл: как data.frame(sdata@data, m2012=r.mean)узнать, какому полигону присвоить какое значение?
Stophface
чтобы правильно добавить r.mean в sdata, мне пришлось сначала удалить список r.mean: r.mean <- unlist (lapply (r.vals, FUN = mean))
cmbarbu
6

Прочитать шейп- файл в SpatialPolygonsDataFrame( readOGRфункция из пакета rgdal)

Считать растр в Rasterобъект ( rasterфункция из пакета raster)

Используйте, extract(raster, spdf)чтобы получить ячейки сетки под каждым полигоном. Тогда беги meanна них.

Повторите с вашим набором растровых изображений ...

Spacedman
источник
как записать шейп-файл (со значениями, извлеченными из растрового изображения)?
Stophface