From thothal at sbox.tugraz.at  Tue Jul  7 14:28:46 2009
From: thothal at sbox.tugraz.at (Thorn Thaler)
Date: Tue, 07 Jul 2009 14:28:46 +0200
Subject: [Rd] Speed up code, profiling, optimization, lapply vs. loops
In-Reply-To: <66f3bd910907070454v52f96af9o1932fe37872d1e20@mail.gmail.com>
References: <4A51B526.7070705@sbox.tugraz.at>	
	<7BDFF633-8E5D-4330-8583-35743292C8F4@stat.berkeley.edu>
	<66f3bd910907070454v52f96af9o1932fe37872d1e20@mail.gmail.com>
Message-ID: <4A533F7E.8010609@sbox.tugraz.at>

Thanks for ur good advices. I've run Rprof and figured out that the 
bottleneck lies within the deviance function, which uses "integrate" to 
determine the deviance numerically. Since I'm not aware of a closed form 
solution for the integral

\int 1/(t^k*(1-t)^l) dt, \forall k,l \in R  ,

I've to rely on the numerical procedure. 'glm' itself is according to 
Rprof suprisingly not that time consuming and I could narrow the lack of 
performance down to multiple calls of the (numerically determined) 
deviance function, which is currently implemented by:

devf <- function(y, mu, ...) {
   mapply(function(y, mu, ...) {
     integrand <- function(u, y, ...)
       return((y - u) / varf(u, ...))
     int <- try(suppressWarnings(
       integrate(integrand, lower=y, upper=mu, y=y, ...)), silent=TRUE)
     if (class(int) == "try-error") {
       stop("Deviance could not be evaluated!\n", int)
     } else {
       return(-2 * int$value)
     }
   }, y, mu, ...)
}

where 'varf' denotes the particular variance function depending on some 
parameter vector.

I'd appreciate any comments on the code. Is it possible to speed up the 
computation?



[Rd] Speed up code, profiling, optimization, lapply vs. loops
Thorn Thaler