Parallelisierte Optimierung in R

Question

Ich renne R auf Linux-Box, die 8 Multicore-Prozessoren hat, und habe ein Optimierungsproblem, das ich möchte durch Parallelisierung der Optimierungsroutine selbst beschleunigen. Wichtig ist, dass dieses Problem beinhaltet (1) mehrere Parameter und (2) inhärent langsam Modell läuft. Ein ziemlich häufiges Problem!

Wer kennt einen parallelisierten Optimierer für solche Gelegenheiten?

Genauer gesagt, führen Löser wie nlm() mehrere Modellauswertungen (zwei pro Parameterwert) jedes Mal durch, wenn der Algorithmus einen Schritt im Parameterraum einlegt, also würde die Parallelisierung dieser Instanz mehrerer Modellläufe die Situation in diesen Situationen erheblich beschleunigen Einige Parameterwerte werden angepasst.

Es scheint, wie Code, der Verwendung des Pakets macht parallel in einer Weise geschrieben werden konnte, dass der Benutzer würde minimal Code-Änderung tun müssen, um aus zu bewegen nlm() oder optim() dieser parallelisierten Optimierungsroutine verwendet wird. Das heißt, man könnte diese Routinen grundsätzlich ohne Änderungen schreiben, außer dass der Schritt des mehrfachen Aufrufs des Modells, wie es bei gradientenbasierten Methoden üblich ist, parallel erfolgen würde.

Idealerweise so etwas wie nlmPara() würde Code nehmen, die

fit <- nlm(MyObjFunc, params0); 

und erfordern nur geringen Modifikationen wie

aussieht, zum Beispiel

fit <- nlmPara(MyObjFunc, params0, ncores=6); 

Gedanken/Anregungen?

PS: Ich habe Schritte unternommen, um diese Modellläufe zu beschleunigen, aber sie sind aus verschiedenen Gründen langsam (d. H. Ich brauche keinen Rat zur Beschleunigung der Modellläufe! ;-)).

Answer 1

Ich verwendete das Paket doSNOW, um einen Code auf 8 Kernen auszuführen. Ich kann einfach & kopieren den Teil des Codes, der auf dieses Paket verweist. Hoffe es hilft!

# use multicore libraries 
     # specify number of cores to use 
    cores<- 8 
     cluster <- makeCluster(cores, type="SOCK") 
     registerDoSNOW(cluster) 

     # check how many cores will be used 
     ncores <- getDoParWorkers() 
    print(paste("Computing algorithm for ", cores, " cores", sep="")) 
     fph <- rep(-100,12) 

     # start multicore cicle on 12 subsets 
     fph <- foreach(i=1:12, .combine='c') %dopar% { 
     PhenoRiceRun(sub=i, mpath=MODIS_LOCAL_DIR, masklocaldir=MASK_LOCAL_DIR, startYear=startYear, tile=tile, evismoothopt=FALSE) 
     } 


    stopCluster(cluster) # check if gives error 
    gc(verbose=FALSE) 

Answer 2

Hier ist eine grobe Lösung, die zumindest etwas verspricht. Vielen Dank an Ben Bolker für den Hinweis, dass viele/die meisten Optimierungsroutinen benutzerspezifische Gradientenfunktionen erlauben.

Ein Testproblem mit mehr Parameterwerten zeigt möglicherweise größere Verbesserungen, aber auf einer 8-Kern-Maschine dauert der Lauf mit der parallelisierten Gradientenfunktion etwa 70% so lange wie die serielle Version. Beachten Sie, dass die hier verwendete grobe Gradientennäherung die Konvergenz zu verlangsamen scheint und somit dem Prozess etwas Zeit hinzufügt.

## Set up the cluster 
require("parallel"); 
.nlocalcores = NULL; # Default to "Cores available - 1" if NULL. 
if(is.null(.nlocalcores)) { .nlocalcores = detectCores() - 1; } 
if(.nlocalcores < 1) { print("Multiple cores unavailable! See code!!"); return()} 
print(paste("Using ",.nlocalcores,"cores for parallelized gradient computation.")) 
.cl=makeCluster(.nlocalcores); 
print(.cl) 


# Now define a gradient function: both in serial and in parallel 
mygr <- function(.params, ...) { 
    dp = cbind(rep(0,length(.params)),diag(.params * 1e-8)); # TINY finite difference 
    Fout = apply(dp,2, function(x) fn(.params + x,...));  # Serial 
    return((Fout[-1]-Fout[1])/diag(dp[,-1]));    # finite difference 
} 

mypgr <- function(.params, ...) { # Now use the cluster 
    dp = cbind(rep(0,length(.params)),diag(.params * 1e-8)); 
    Fout = parCapply(.cl, dp, function(x) fn(.params + x,...)); # Parallel 
    return((Fout[-1]-Fout[1])/diag(dp[,-1]));     # 
} 


## Lets try it out! 
fr <- function(x, slow=FALSE) { ## Rosenbrock Banana function from optim() documentation. 
    if(slow) { Sys.sleep(0.1); } ## Modified to be a little slow, if needed. 
    x1 <- x[1] 
    x2 <- x[2] 
    100 * (x2 - x1 * x1)^2 + (1 - x1)^2 
} 

grr <- function(x, slow=FALSE) { ## Gradient of 'fr' 
    if(slow) { Sys.sleep(0.1); } ## Modified to be a little slow, if needed. 
    x1 <- x[1] 
    x2 <- x[2] 
    c(-400 * x1 * (x2 - x1 * x1) - 2 * (1 - x1), 
    200 *  (x2 - x1 * x1)) 
} 

## Make sure the nodes can see these functions & other objects as called by the optimizer 
fn <- fr; # A bit of a hack 
clusterExport(cl, "fn"); 

# First, test our gradient approximation function mypgr 
print(mypgr(c(-1.2,1)) - grr(c(-1.2,1))) 

## Some test calls, following the examples in the optim() documentation 
tic = Sys.time(); 
fit1 = optim(c(-1.2,1), fr, slow=FALSE);       toc1=Sys.time()-tic 
fit2 = optim(c(-1.2,1), fr, gr=grr, slow=FALSE, method="BFGS"); toc2=Sys.time()-tic-toc1 
fit3 = optim(c(-1.2,1), fr, gr=mygr, slow=FALSE, method="BFGS"); toc3=Sys.time()-tic-toc1-toc2 
fit4 = optim(c(-1.2,1), fr, gr=mypgr, slow=FALSE, method="BFGS"); toc4=Sys.time()-tic-toc1-toc2-toc3 


## Now slow it down a bit 
tic = Sys.time(); 
fit5 = optim(c(-1.2,1), fr, slow=TRUE);       toc5=Sys.time()-tic 
fit6 = optim(c(-1.2,1), fr, gr=grr, slow=TRUE, method="BFGS"); toc6=Sys.time()-tic-toc5 
fit7 = optim(c(-1.2,1), fr, gr=mygr, slow=TRUE, method="BFGS"); toc7=Sys.time()-tic-toc5-toc6 
fit8 = optim(c(-1.2,1), fr, gr=mypgr, slow=TRUE, method="BFGS"); toc8=Sys.time()-tic-toc5-toc6-toc7 

print(cbind(fast=c(default=toc1,exact.gr=toc2,serial.gr=toc3,parallel.gr=toc4), 
      slow=c(toc5,toc6,toc7,toc8))) 

Answer 3

Wie Sie eine Antwort nicht akzeptiert haben, könnte diese Idee helfen: Für globale Optimierung des Paket DEoptim() eine eingebaute Option für die parallele Optimierung hat. Gute Sache ist, es ist einfach zu bedienen und Dokumentation gut geschrieben.

c.f. http://www.jstatsoft.org/v40/i06/paper (derzeit unten)

http://cran.r-project.org/web/packages/DEoptim/index.html

Vorsicht: Differenz Evolglobal Optimizern noch in Einheimischen laufen könnten.