Обработка вектора в растре быстрее с помощью R

9

Я конвертирую вектор в растр в R. Однако процесс был слишком долгим. Есть ли возможность поместить скрипт в многопоточную или графическую обработку, чтобы сделать это быстрее?

Мой скрипт в растеризованный вектор.

r.raster = raster()
extent(r.raster) = extent(setor) #definindo o extent do raster
res(r.raster) = 10 #definindo o tamanho do pixel
setor.r = rasterize(setor, r.raster, 'dens_imov')

r.raster

класс: RasterLayer размеры: 9636, 11476, 110582736 (nrow, ncol, ncell) разрешение: 10, 10 (x, y) экстент: 505755, 620515, 8555432, 8651792 (координаты xmin, xmax, ymin, ymax). ссылка : + proj = longlat + datum = WGS84 + ellps = WGS84 + towgs84 = 0,0,0

Setor

Класс: SpatialPolygonsDataFrame Особенности: 5419 экстент: 505755, 620515.4, 8555429, 8651792 (xmin, xmax, ymin, ymax) координат. ссылка : + proj = utm + zone = 24 + south + ellps = GRS80 + unit = m + no_defs переменные: 6 имен: ID, CD_GEOCODI, TIPO, dens_imov, area_m, domicilios1 min значения: 35464, 290110605000001, RURAL, 0,00000003,100004, Максимальные значения 1,0000: 58468, 293320820000042, URBANO, 0,54581673,99996, 99,0000

Печать сетора введите описание изображения здесь

Диого Карибе
источник
Можете ли вы опубликовать резюме сетора и r.raster? Я хотел бы иметь некоторое представление о количестве объектов в setor и размерах r.raster. просто напечатайте их хорошо
mdsumner
Я поставил резюме в теле вопроса.
Диого Карибе
Не резюме, просто распечатай - информация, которую я попросил нас, не tgere
mdsumner
Извините, я поставил печать.
Диого Карибе
Ах, разочарован, я не думал об этом до тех пор, пока не увидел распечатку - убедитесь, что проекция растра совпадает с полигонами, в данный момент это не так - попробуйте r <- растр (setor); res (r) <- 10; setor.r = rasterize (setor, r, 'dens_imov') - но также попробуйте сначала установить res (r) <- 250, чтобы вы
знали,

Ответы:

17

Я попытался «распараллелить» функцию, rasterizeиспользуя Rпакет parallelследующим образом:

  1. разбить объект SpatialPolygonsDataFrame на nчасти
  2. rasterize каждая часть отдельно
  3. объединить все части в один растр

В моем компьютере rasterizeфункция распараллеливания выполнялась в 2,75 раза меньше, чем rasterizeфункция без распараллеливания .

Примечание. Приведенный ниже код загружает шейп-файл многоугольника (~ 26,2 МБ) из Интернета. Вы можете использовать любой объект SpatialPolygonDataFrame. Это только пример.

Загрузите библиотеки и пример данных:

# Load libraries
library('raster')
library('rgdal')

# Load a SpatialPolygonsDataFrame example
# Load Brazil administrative level 2 shapefile
BRA_adm2 <- raster::getData(country = "BRA", level = 2)

# Convert NAMES level 2 to factor 
BRA_adm2$NAME_2 <- as.factor(BRA_adm2$NAME_2)

# Plot BRA_adm2
plot(BRA_adm2)
box()

# Define RasterLayer object
r.raster <- raster()

# Define raster extent
extent(r.raster) <- extent(BRA_adm2)

# Define pixel size
res(r.raster) <- 0.1

BrazilSPDF

Рисунок 1: График Бразилии SpatialPolygonsDataFrame

Пример простой темы

# Simple thread -----------------------------------------------------------

# Rasterize
system.time(BRA_adm2.r <- rasterize(BRA_adm2, r.raster, 'NAME_2'))

Время в моем ноутбуке:

# Output:
# user  system elapsed 
# 23.883    0.010   23.891

Пример многопоточности

# Multithread -------------------------------------------------------------

# Load 'parallel' package for support Parallel computation in R
library('parallel')

# Calculate the number of cores
no_cores <- detectCores() - 1

# Number of polygons features in SPDF
features <- 1:nrow(BRA_adm2[,])

# Split features in n parts
n <- 50
parts <- split(features, cut(features, n))

# Initiate cluster (after loading all the necessary object to R environment: BRA_adm2, parts, r.raster, n)
cl <- makeCluster(no_cores, type = "FORK")
print(cl)

# Parallelize rasterize function
system.time(rParts <- parLapply(cl = cl, X = 1:n, fun = function(x) rasterize(BRA_adm2[parts[[x]],], r.raster, 'NAME_2')))

# Finish
stopCluster(cl)

# Merge all raster parts
rMerge <- do.call(merge, rParts)

# Plot raster
plot(rMerge)

BrazilRaster

Рисунок 2: График растров Бразилии

Время в моем ноутбуке:

# Output:
# user  system elapsed 
# 0.203   0.033   8.688 

Больше информации о распараллеливании в R :

Гусман
источник
Очень хороший ответ!
Диого Карибе
Разве вы не просто задали число ядер на машине?
Сэм
@ Сам, я думаю, это должно работать без проблем, но я не знаю, лучше это или нет! Я предположил, что если бы я разделил функции на n частей, равных количеству ядер, возможно, одну из этих частей можно было бы проще обрабатывать, а ядро, которое обрабатывало его, было бы бесполезным! Однако, если у вас больше деталей, чем ядер, когда одно ядро ​​завершит обработку одной части, это займет другую часть. Но, конечно, я не уверен! Любая помощь по этому вопросу будет принята с благодарностью.
Гусман
Я собираюсь провести несколько тестов сегодня вечером. На небольшом шейп-файле (примерно 25 на 25 км), растеризованном до 50 м, есть небольшое улучшение в использовании n = 2,4 или 8 против n = 20, 30 или до 50. Сегодня вечером я добавлю очень большой шейп-файл и растеризировать до 25м. Обработка одного ядра занимает 10 часов, поэтому мы посмотрим, что делают различные значения n !! (n = 50 - чуть меньше 1 часа)
Сэм
@ Guzmán Я снова запускаю код. Тем не менее, это исправило некоторые ошибки и не знаю почему. Вы можете помочь мне? Ошибка в checkForRemoteErrors (val): 7 узлов выдали ошибки; первая ошибка: объект 'BRA_adm2' не найден
Diogo Caribé