Проблема

Я пытаюсь использовать модель, которая требует от меня преобразования очень больших растров (около 10 миллионов ячеек) в data.frame, я делаю это на общем сервере, и поэтому я пытаюсь оптимизировать, чтобы уменьшить Используемая оперативная память и, надеюсь, не сильно снизит скорость. Для этого я написал 2 функции, но безуспешно. RPUBS с результатами моих попыток находится по этой ссылке, а github с rmd для этого здесь, включая файлы rds для профилирования profvis.

Первая функция

Сначала давайте загрузим пакеты, которые нам понадобятся:

# For spaital manipulation
library(raster)
# For benchmarking speed and memory
library(profvis)
# To parallelize operations
library(doParallel)
# To parallelize operations
library(foreach)
# For combination and looping
library(purrr)
# Data wranggling
library(dplyr)
library(data.table)

плитка

Основной способ уменьшить использование оперативной памяти — вместо обработки одного большого растра преобразовать его в мозаичные растры, для чего я разработал следующую функцию:

SplitRas <- function(Raster,ppside, nclus = 1){
  TempRast <- paste0(getwd(), "/Temp")
  h        <- ceiling(ncol(Raster)/ppside)
  v        <- ceiling(nrow(Raster)/ppside)
  agg      <- aggregate(Raster,fact=c(h,v))
  agg[]    <- 1:ncell(agg)
  agg_poly <- rasterToPolygons(agg)
  names(agg_poly) <- "polis"
  r_list <- list()
  if(nclus == 1){
    for(i in 1:ncell(agg)){
      dir.create(TempRast)
      rasterOptions(tmpdir=TempRast)
      e1          <- extent(agg_poly[agg_poly$polis==i,])
      r_list[[i]] <- crop(Raster,e1)
      if((freq(r_list[[i]], value=NA)/ncell(r_list[[i]])) != 1){
        writeRaster(r_list[[i]],filename=paste("SplitRas",i,sep=""),
                  format="GTiff",datatype="FLT4S",overwrite=TRUE)
      }
      unlink(TempRast, recursive = T, force = T)
    }
  } else if(nclus != 1){
    cl <- parallel::makeCluster(nclus)
    doParallel::registerDoParallel(cl)
    r_list <- foreach(i = 1:ncell(agg), .packages = c("raster")) %dopar% {
      dir.create(TempRast)
      rasterOptions(tmpdir=TempRast)
      e1          <- extent(agg_poly[agg_poly$polis==i,])
      Temp <- crop(Raster,e1)
      if((raster::freq(Temp, value=NA)/ncell(Temp)) != 1){
        writeRaster(Temp,filename=paste("SplitRas",i,sep=""),
                    format="GTiff",datatype="FLT4S",overwrite=TRUE)
      }
      unlink(TempRast, recursive = T, force = T)
      Temp
    }
    parallel::stopCluster(cl)
  }
}

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

Эта функция имеет 3 аргумента:

• Растр. Стопка, которую вы хотите разделить на плитки.
• ppside: количество горизонтальных и вертикальных плиток, которые вы получите в итоге, окончательное количество плиток будет ppside*ppside, поэтому, если ppside равно 3, у вас будет 9 плиток.
• nclus: количество кластеров, которое вы будете использовать для распараллеливания и ускорения обработки.

В конце этой функции у вас будет ppside*ppside количество плиток, сохраненных в вашем рабочем каталоге, и все они будут называться SplitRasN.tif, где N — номер плитки. В качестве примера воспользуемся биоклиматическими переменными, доступными в растровом пакете:

Bios <- getData('worldclim', var='bio', res=10)

Я протестирую это преобразование в другое количество плиток, как показано на следующем рисунке.

Преобразование растра в тайлы, а затем из тайлов во фрейм данных (пример)

поэтому сначала мы будем использовать функцию SplitRas, чтобы получить 16 плиток, используя 4 ядра:

SplitRas(Raster = Bios, ppside = 4, nclus = 4)

Это вернет следующие файлы: r list.files(pattern = "SplitRas")

Чтобы собрать эти плитки в один фрейм данных со всеми ячейками, не относящимися к NA, нам нужен список плиток, который мы получаем с помощью:

Files <- list.files(pattern = "SplitRas", full.names = T)

Который мы можем использовать затем в следующей функции:

SplitsToDataFrame <- function(Splits, ncores = 1){
  TempRast <- paste0(getwd(), "/Temp")
  if(ncores == 1){
    Temps <- list()
    for(i in 1:length(Splits)){
      dir.create(TempRast)
      rasterOptions(tmpdir=TempRast)
      Temp <- stack(Splits[i])
      Temp <- as.data.frame(Temp, row.names = NULL, col.names = NULL, xy =TRUE)
      colnames(Temp)[3:ncol(Temp)] <- paste0("Var", 1:length(3:ncol(Temp)))
      Temps[[i]] <- Temp[complete.cases(Temp),]
      gc()
      unlink(TempRast, recursive = T, force = T)
      message(i)
    }
    Temps <- data.table::rbindlist(Temps)
  } else if(ncores > 1){
    cl <- parallel::makeCluster(ncores)
    doParallel::registerDoParallel(cl)
    Temps <- foreach(i = 1:length(Splits), .packages = c("raster", "data.table")) %dopar%{
      dir.create(TempRast)
      rasterOptions(tmpdir=TempRast)
      Temp <- stack(Splits[i])
      Temp <- as.data.frame(Temp, row.names = NULL, col.names = NULL, xy =TRUE)
      colnames(Temp)[3:ncol(Temp)] <- paste0("Var", 1:length(3:ncol(Temp)))
      gc()
      unlink(TempRast, recursive = T, force = T)
      Temp[complete.cases(Temp),]
    }
    Temps <- data.table::rbindlist(Temps)
    parallel::stopCluster(cl)
  }
  return(Temps)
}

Где Splits — это векторы символов с путями к тайлам, а ncores — количество ядер, используемых для распараллеливания. Это приведет к кадру данных с непустыми ячейками.

DF <- SplitsToDataFrame(Splits = Files, ncores = 4)

Тестирование памяти (что я пробовал)

Сначала я сгенерировал плитки для 1, 2, 4, 8, 10 и 12 сторон.

sides <- c(1,2,4,8,10, 12)

Home <- getwd()
AllFiles <- list()
for(i in 1:length(sides)){
  dir.create(path = paste0(Home, "/Sides_", sides[i]))
  setwd(paste0(Home, "/Sides_", sides[i]))
  SplitRas(Raster = Bios, ppside = sides[i], nclus = ifelse(sides[i] < 4, sides[i], 4))
  AllFiles[[i]] <- list.files(pattern = "SplitRas", full.names = T) %>% stringr::str_replace_all("\\./", paste0(getwd(), "/"))
}

setwd(Home)

А затем профилировал функцию, используя только последовательный вариант функции

library(profvis)
P <- profvis({
    P1 <- SplitsToDataFrame(Splits = AllFiles[[1]])
    P2 <- SplitsToDataFrame(Splits = AllFiles[[2]])
    P3 <- SplitsToDataFrame(Splits = AllFiles[[3]])
    P4 <- SplitsToDataFrame(Splits = AllFiles[[4]])
    P5 <- SplitsToDataFrame(Splits = AllFiles[[5]])
})
P
htmlwidgets::saveWidget(P, "profile.html")
saveRDS(P, "P.rds")

Однако, как видно на следующем рисунке (более подробную информацию можно найти в ссылках Rpubs выше), оперативная память более или менее такая же, как и раньше, но время увеличилось (эта последняя часть ожидалась).

.

Некоторые дополнительные вещи

Наконец, когда я попытался сравнить использование ОЗУ параллельно, мне показалось, что profvis этого не фиксирует. Любая идея о том, как это проверить?

library(profvis)
PPar <- profvis({
    P1 <- SplitsToDataFrame(Splits = AllFiles[[3]], ncores = 1)
    P2 <- SplitsToDataFrame(Splits = AllFiles[[3]], ncores = 2)
    P3 <- SplitsToDataFrame(Splits = AllFiles[[3]], ncores = 4)
    P4 <- SplitsToDataFrame(Splits = AllFiles[[3]], ncores = 7)
})
PPar
htmlwidgets::saveWidget(PPar, "profileParallel.html")
saveRDS(PPar, "PPar.rds")