How about? > fun <- function(u, m, s) 1/ (1 + exp(- (B[1] + B[2] * + (m + u)))) * dnorm(u, 0, s) > eta <- sapply(1:nrow(X), function(i) integrate(fun, + -Inf, Inf, m=x[i], s=sem1[i])$value) > eta [1] 0.3018974 0.6780813 0.4804123 0.4845773 0.3886028 0.3336432 0.2541276 [8] 0.6894951 0.6649838 0.2385417
---------------------------------------------- David L Carlson Associate Professor of Anthropology Texas A&M University College Station, TX 77843-4352 > -----Original Message----- > From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- > project.org] On Behalf Of Doran, Harold > Sent: Friday, December 07, 2012 9:14 AM > To: David Winsemius > Cc: r-help@r-project.org > Subject: Re: [R] Vectorizing integrate() > > David et al > > Thanks, I should have made the post more complete. I routinely use > apply functions, but often avoid mapply() as I find it so non- > intuitive. In this instance, I think the situation demands I change > that position, though. > > Reproducible code for the current implementation of the function is > > B <- c(0,1) > sem1 = runif(10, 1, 2) > x <- rnorm(10) > X <- cbind(1, x) > eta <- numeric(10) > > for(j in 1:nrow(X)){ > fun <- function(u) 1/ (1 + exp(- (B[1] + B[2] * (x[j] + u)))) * > dnorm(u, 0, sem1[j]) > eta[j] <- integrate(fun, -Inf, Inf)$value > } > > I can't get my head around how mapply() would work here. It accepts as > its first argument a function. But, in my case I have two functions: > the user defined integrand, fun(), an then of course calling the R > function integrate(). > > I was thinking maybe along these lines, but this is obviously wrong. > > mapply(integrate(function(u) 1/ (1 + exp(- (B[1] + B[2] * (x + u)))) * > dnorm(u, 0, sem1), -Inf, Inf)$value, MoreArgs = list(B, x, sem1)) > > > > > -----Original Message----- > > From: David Winsemius [mailto:dwinsem...@comcast.net] > > Sent: Thursday, December 06, 2012 1:59 PM > > To: Doran, Harold > > Cc: r-help@r-project.org > > Subject: Re: [R] Vectorizing integrate() > > > > > > On Dec 6, 2012, at 10:10 AM, Doran, Harold wrote: > > > > > I have written a program to solve a particular logistic regression > problem > > using IRLS. In one step, I need to integrate something out of the > linear > > predictor. The way I'm doing it now is within a loop and it is as you > would > > expect slow to process, especially inside an iterative algorithm. > > > > > > I'm hoping there is a way this can be vectorized, but I have not > found > > > it so far. The portion of code I'd like to vectorize is this > > > > > > for(j in 1:nrow(X)){ > > > fun <- function(u) 1/ (1 + exp(- (B[1] + B[2] * (x[j] + u)))) * > dnorm(u, 0, > > sd[j]) > > > eta[j] <- integrate(fun, -Inf, Inf)$value } > > > > > > > The Vectorize function is just a wrapper to mapply. If yoou are able > to get > > that code to execute properly for your un-posted test cases, then why > not > > use mapply? > > > > > > > Here X is an n x p model matrix for the fixed effects, B is a > vector with the > > estimates of the fixed effects at iteration t, x is a predictor > variable in the jth > > row of X, and sd is a variable corresponding to x[j]. > > > > > > Is there a way this can be done without looping over the rows of X? > > > > > > Thanks, > > > Harold > > > > > > [[alternative HTML version deleted]] > > > > > > ______________________________________________ > > > R-help@r-project.org mailing list > > > https://stat.ethz.ch/mailman/listinfo/r-help > > > PLEASE do read the posting guide > > > http://www.R-project.org/posting-guide.html > > > and provide commented, minimal, self-contained, reproducible code. > > > > David Winsemius, MD > > Alameda, CA, USA > > ______________________________________________ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting- > guide.html > and provide commented, minimal, self-contained, reproducible code. ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.