Memproses vektor ke raster lebih cepat dengan R

9

Saya mengubah vektor menjadi raster di R. Namun prosesnya terlalu lama. Apakah ada kemungkinan untuk menempatkan skrip ke dalam pemrosesan multithread atau GPU untuk melakukannya lebih cepat?

Script saya untuk vektor raster.

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')

raster

kelas: Dimensi RasterLayer: 9636, 11476, 110582736 (nrow, ncol, ncell) resolusi: 10, 10 (x, y) sejauh: 505755, 620515, 8555432, 8651792 (xmin, xmax, ymin, ymax) coord. ref. : + proj = longlat + datum = WGS84 + ellps = WGS84 + towgs84 = 0,0,0

setor

kelas: SpatialPolygonsDataFrame fitur: 5419 luas: 505755, 620515.4, 8555429, 8651792 (xmin, xmax, ymin, ymax) coord. ref. : + proj = utm + zona = 24 + selatan + ellps = GRS80 + unit = m + variabel no_defs: 6 nama: ID, CD_GEOCODI, TIPO, dens_imov, area_m, domicilios1 nilai minimum: 35464, 290110605000001, RURAL, 0,00000003,100004, Nilai maks 1,0000: 58468, 293320820000042, URBANO, 0,54581673,99996, 99,0000

Cetak setor masukkan deskripsi gambar di sini

raster r rasterization rgdal parallel Diogo Caribé
sumber
Bisakah Anda memposting ringkasan setor dan r.raster? Saya ingin memiliki beberapa gagasan tentang jumlah objek dalam setor dan dimensi r.raster. cukup cetak saja
mdsumner
Saya menaruh ringkasan di tubuh pertanyaan.
Diogo Caribé
Bukan ringkasan, cukup cetak - info yang saya minta tidak kami tgere
mdsumner
Maaf, saya pasang cetakannya.
Diogo Caribé
Ah, kecewa saya tidak memikirkan hal ini sampai saya melihat cetakannya - pastikan proyeksi raster cocok dengan poligon, saat ini tidak - coba r <- raster (setor); res (r) <- 10; setor.r = rasterize (setor, r, 'dens_imov') - tetapi juga coba, atur res (r) <- 250 terlebih dahulu sehingga Anda mendapatkan gagasan tentang berapa lama versi resolusi tinggi akan memakan waktu
mdsumner

Jawaban:

17

Saya mencoba "memparalelkan" fungsi rasterizemenggunakan Rpaket paralleldengan cara ini:

  1. pisahkan objek SpatialPolygonsDataFrame menjadi beberapan bagian
  2. rasterize setiap bagian secara terpisah
  3. menggabungkan semua bagian menjadi satu raster

Di komputer saya, rasterizefungsi yang diparalelkan mengambil 2,75 kali lebih sedikit dari fungsi yang tidak diparalelkan rasterize.

Catatan: kode di bawah ini mengunduh shapefile poligon (~ 26,2 MB) dari web. Anda bisa menggunakan objek SpatialPolygonDataFrame apa pun. Ini hanya sebuah contoh.

Muat perpustakaan dan contoh data:

# 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

Gambar 1: Plot Brasil SpatialPolygonsDataFrame

Contoh utas sederhana

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

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

Waktu di laptop saya:

# Output:
# user  system elapsed 
# 23.883    0.010   23.891

Contoh thread multithread

# 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

Gambar 2: Plot Raster Brasil

Waktu di laptop saya:

# Output:
# user  system elapsed 
# 0.203   0.033   8.688

Info lebih lanjut tentang paralelisasi di R :

Guzmán
sumber
Jawaban yang sangat bagus
Diogo Caribé
Apakah Anda tidak menetapkan n sebagai jumlah inti pada mesin?
Sam
@ Sam Saya pikir itu harus bekerja tanpa masalah tapi saya tidak tahu apakah ini lebih baik atau tidak! Saya berasumsi bahwa jika saya memisahkan fitur dalam n bagian yang sama dengan jumlah core, mungkin salah satu dari bagian ini dapat lebih mudah diproses dan inti yang memprosesnya akan tanpa digunakan! Namun, jika Anda memiliki lebih banyak bagian dari inti ketika satu inti selesai memproses satu bagian, itu akan mengambil bagian lainnya. Tapi yang pasti, saya tidak yakin! Bantuan apa pun tentang masalah ini akan dihargai.
Guzmán
Saya akan menjalankan beberapa tes malam ini. Pada shapefile kecil (kira-kira 25 km x 25 km), dirasterisasi menjadi 50 m, ada sedikit peningkatan dalam menggunakan n = 2,4 atau 8 melawan n = 20, 30 atau lebih dari 50. Saya akan masuk dalam shapefile yang sangat besar malam ini dan rasterize menjadi 25m. Pemrosesan inti tunggal adalah 10 jam sehingga kita akan melihat apa nilai n yang berbeda lakukan !! (n = 50 tepat di bawah 1 jam)
Sam
@ Guzmán Saya menjalankan kode lagi. Namun, beberapa kesalahan muncul kembali dan tidak tahu mengapa. Bisakah kamu membantuku? Kesalahan dalam checkForRemoteErrors (val): 7 node menghasilkan kesalahan; kesalahan pertama: objek 'BRA_adm2' tidak ditemukan
Diogo Caribé