effiziente Nutzung von Rasterfunktionen in r

Question

Ich habe 500 + Punkte in einem SpatialPointsDataFrame Objekt; Ich habe ein 1,7 GB (200.000 Zeilen x 200.000 Spalten) raster Objekt. Ich möchte eine Tabelle der Werte der Rasterzellen innerhalb eines Puffers um jeden der 500+ Punkte haben.

Ich habe es geschafft, mit dem Code unten zu erreichen (ich habe eine Menge Inspiration bekommen from here.). Es ist jedoch langsam zu laufen und ich möchte es schneller laufen lassen. Bei Puffern mit "kleinen" Breiten, z. B. 5 km/h, sogar 15 km (~ 1 Mio. Zellen), läuft es sogar OK, aber es wird sehr langsam, wenn der Puffer auf 100 km ansteigt (~ 42 Mio. Zellen).

Ich könnte leicht auf die Schleife unten verbessern, indem Sie etwas aus der apply Familie und/oder eine parallele Schleife verwenden. Aber mein Verdacht ist, dass es langsam ist, weil das raster Paket 400Mb + temporäre Dateien für jede Interaktion der Schleife schreibt.

# packages 
library(rgeos) 
library(raster) 
library(rgdal) 

myPoints = readOGR(points_path, 'myLayer') 
myRaster = raster(raster_path) 

myFunction = function(polygon_obj, raster_obj) { 
    # this function return a tabulation of the values of raster cells 
    # inside a polygon (buffer) 

    # crop to extent of polygon 
    clip1 = crop(raster_obj, extent(polygon_obj)) 

    # crops to polygon edge & converts to raster 
    clip2 = rasterize(polygon_obj, clip1, mask = TRUE) 

    # much faster than extract 
    ext = getValues(clip2) 

    # tabulates the values of the raster in the polygon 
    tab = table(ext) 

    return(tab) 
} 

# loop over the points 
ids = unique(myPoints$ID) 
for (id in ids) { 

    # select point 
    myPoint = myPoints[myPoints$ID == id, ] 

    # create buffer 
    myPolygon = gBuffer(spgeom = myPoint, byid = FALSE, width = myWidth) 

    # extract the data I want (projections, etc are fine) 
    tab = myFunction(myPolygon, myRaster) 

    # do stuff with tab ... 
} 

Meine Fragen:

1. Habe ich recht teilweise die Schreiboperationen schuld? Wenn ich all diese Schreiboperationen vermeiden könnte, würde dieser Code schneller laufen? Ich habe Zugriff auf eine Maschine mit 32 GB RAM - also kann ich davon ausgehen, dass ich die raster in den Speicher laden kann und keine temporären Dateien schreiben muss?

2. Was könnte ich noch tun, um die Effizienz in diesem Code zu verbessern?

Answer 1

Ich glaube, Sie es wie dieser

library(raster) 
library(rgdal) 
myPoints <- readOGR(points_path, 'myLayer') 
myRaster <- raster(raster_path) 
e <- extract(myRaster, myPoints, buffer=myWidth) 

Und dann so etwas wie

etab <- sapply(e, table) 

Es ist schwer nähern sollte Ihre Frage # 1, wie wir wissen, nicht genug zu beantworten über Ihre Daten (wir wissen nicht, wie viele Zellen durch einen Puffer "100 km" abgedeckt sind). Sie können jedoch festlegen, wann in die Datei mit der Funktion rasterOptions geschrieben werden soll. Sie bemerken, dass getValues schneller ist als Auszug, basierend auf dem Beitrag, den Sie verlinken, aber ich denke, das ist falsch oder zumindest nicht sehr wichtig. Die Kombination von crop, und getValues sollte eine ähnliche Leistung wie extract (die fast genau das unter der Haube tut) haben. Wenn Sie diese Route dennoch verwenden, sollten Sie eine leere RasterLayer übergeben, die von raster(myRaster) erstellt wird, um schneller zu beschneiden.