提高性能/速度

Question

我需要從1303柵格（每個柵格有1個月的數據）中獲取數據，併爲柵格中的每個柵格單元創建一個時間序列。最後，我會將所有時間序列合併成一個大型（動物園）文件。

我有代碼可以做到這一點（我嘗試了一小部分數據集，它的工作原理），但它似乎只是堆疊柵格（現在超過2小時，仍然在計數）和這不是較慢的部分，這將是時間序列。所以這裏是我的代碼，如果有人知道更快的方式來堆疊柵格和/或創建時間序列（也許沒有雙循環？）請幫助...

我不知道任何其他編程語言，但這對R來說太過分了嗎？

files <- list.files(pattern=".asc") 
pat <- "^.*pet_([0-9]{1,})_([0-9]{1,}).asc$" 
ord_files <- as.Date(gsub(pat, sprintf("%s-%s-01", "\\1", "\\2"), files)) 
files<-files[order(ord_files)] 


#using "raster" package to import data 
s<- raster(files[1]) 
pet<-vector() 
for (i in 2:length(files)) 
{ 
r<- raster(files[i]) 
s <- stack(s, r) 
} 

#creating a data vector 
beginning = as.Date("1901-01-01") 
full <- seq(beginning, by='1 month', length=length(files)) 
dat<-as.yearmon(full) 

#building the time series 
for (lat in 1:360) 
for (long in 1:720) 
{ 
pet<-as.vector(s[lat,long]) 
x <- xts(pet, dat) 
write.zoo(x,file=paste("P:/WRSRL/Users1/ncgk/IBERIA/cru_pet/zoo/","lat",lat,"long",long,".csv", sep="") , sep=",") 
} 

Answer 1

我將在這裏重新發布我的意見，並給一個更好的例子：

總體思路：在執行「raster'環之前分配對於s的空間。如果將s和r連接到循環中的新對象s，則R必須爲每次迭代爲s分配新的內存。這真的很慢，特別是如果s很大。

s <- c() 
system.time(for(i in 1:1000){ s <- c(s, rnorm(100))}) 
# user system elapsed 
# 0.584 0.244 0.885 

s <- rep(NA, 1000*100) 
system.time(for(i in seq(1,1000*100,100)){ s[i:(i+99)] <- rnorm(100) }) 
# user system elapsed 
# 0.052 0.000 0.050 

正如您所看到的，預分配快了大約10倍。

不幸的是我不熟悉raster和stack，所以我不能告訴你如何將它應用到你的代碼中。

Answer 2

像這樣的東西應該工作（如果你有足夠的內存）：

#using "raster" package to import data 
rlist <- lapply(files, raster) 
s <- do.call(stack, rlist) 
rlist <- NULL # to allow freeing of memory 

它加載所有raster物體進入大名單，然後調用stack一次。

這裏的速度上漲的一個例子：1.25秒VS 8秒60頁的文件 - 但你的舊代碼是二次在時間上的收益更多的文件要高得多......

library(raster) 
f <- system.file("external/test.grd", package="raster") 
files <- rep(f, 60) 

system.time({ 
rlist <- lapply(files, raster) 
s <- do.call(stack, rlist) 
rlist <- NULL # to allow freeing of memory 
}) # 1.25 secs 

system.time({ 
s<- raster(files[1]) 
for (i in 2:length(files)) { 
    r<- raster(files[i]) 
    s <- stack(s, r) 
} 
}) # 8 secs 

Answer 3

第一位可能僅僅是：

s <- stack(files) 

之所以創建堆棧是有點慢是每個文件需要被打開和檢查，看它是否具有相同的nrow，NcoI位等的其他文件。如果你是絕對肯定是這樣的話，你可以使用這樣的快捷方式（一般不推薦）

quickStack <- function(f) { 
r <- raster(f[1]) 
ln <- extension(basename(f), '') 
s <- stack(r) 
s@layers <- sapply(1:length(f), function(x){ r@file@name = f[x]; r@layernames=ln[x]; r@data@haveminmax=FALSE ; r }) 
s@layernames <- ln 
s 
} 

quickStack(files) 

你或許可以也加快第二部分是在下面的例子中，取決於有多少RAM你有。

讀一行一行：

for (lat in 1:360) { 
pet <- getValues(s, lat, 1) 
for (long in 1:720) { 
    x <- xts(pet[long,], dat) 
    write.zoo(x,file=paste("P:/WRSRL/Users1/ncgk/IBERIA/cru_pet/zoo/","lat",lat,"long",long,".csv", sep="") , sep=",") 
} 
} 

更加極端，閱讀一氣呵成的所有值：

pet <- getValues(s) 
for (lat in 1:360) { 
for (long in 1:720) { 
    cell <- (lat-1) * 720 + long 
    x <- xts(pet[cell,], dat) 
    write.zoo(x,file=paste("P:/WRSRL/Users1/ncgk/IBERIA/cru_pet/zoo/","lat",lat,"long",long,".csv", sep="") , sep=",") 
} 
} 

Answer 4

我嘗試另一種方式來處理大量文件。 首先，我將時間序列柵格合併到NetCDF格式的一個文件中， 使用write.Raster（x，format =「CDF」，..） 然後每年只讀一個文件，這次我使用了磚塊（netcdffile，varname =''），它的讀取節省了很多。 但是，我需要根據一些預定義的格式保存所有年份的每個單元格的值，其中我使用write.fwf（x = v，...，append = TRUE） ，但需要很長時間才能使用近500,000點。 任何人都有相同的經驗和幫助如何加快這一進程？ 這裏是我的每個點提取所有的值代碼：

weather4Point <- function(startyear,endyear) 
{ 

    for (year in startyear:endyear) 
    { 
    #get the combined netCDF file 

    tminfile <- paste("tmin","_",year,".nc",sep='') 

    b_tmin <- brick(tminfile,varname='tmin') 

    pptfile <- paste("ppt","_",year,".nc",sep='') 

    b_ppt <- brick(pptfile,varname='ppt') 

    tmaxfile <- paste("tmax","_",year,".nc",sep='') 

    b_tmax <- brick(tmaxfile,varname='tmax') 

    #Get the first year here!!! 

    print(paste("processing year :",year,sep='')) 

    for(l in 1:length(pl)) 
    { 
     v <- NULL 

     #generate file with the name convention with t_n(latitude)w(longitude).txt, 5 digits after point should be work 

     filename <- paste("c:/PRISM/MD/N",round(coordinates(pl[l,])[2],5),"W",abs(round(coordinates(pl[l,])[1],5)),".wth",sep='') 

     print(paste("processing file :",filename,sep=''))    

     tmin <- as.numeric(round(extract(b_tmin,coordinates(pl[l,])),digits=1)) 

     tmax <- as.numeric(round(extract(b_tmax,coordinates(pl[l,])),digits=1)) 

     ppt <- as.numeric(round(extract(b_ppt,coordinates(pl[l,])),digits=2)) 

     v <- cbind(tmax,tmin,ppt) 

     tablename <- c("tmin","tmax","ppt") 

     v <- data.frame(v) 

     colnames(v) <- tablename 

     v["default"] <- 0 

     v["year"] <- year 

     date <- seq(as.Date(paste(year,"/1/1",sep='')),as.Date(paste(year,"/12/31",sep='')),"days") 

     month <- as.numeric(substr(date,6,7)) 

     day <- as.numeric(substr(date,9,10)) 

     v["month"] <- month 

     v["day"] <- day 

     v <- v[c("year","month","day","default","tmin","tmax","ppt")] 

     #write into a file with format 
     write.fwf(x=v,filename,append=TRUE,na="NA",rownames=FALSE,colnames=FALSE,width=c(6,3,3,5,5,5,6)) 
    } 
    } 
}