[R] vectorizing a function

[Robin Hankin]

Hi Bill


> > besselI() does The Right Thing (tm), but it is internal ; what is the
> > best way to vectorize this type of function?
> > 

> > [my clunky hypergeometric function deleted]


> [Bill's disclaimer]

 
> 	(NOTE: this is also vectorized wrt A, B and C as a kind of bonus.)
> 
> 	hyper <- function(A, B, C, z, tol = sqrt(.Machine$double.eps), maxit
> = 10000) {
> 	    term <- z
> 	    term[ ] <- 1
> 	    partial.sum <- term
> 	    n <- 0
> 	    while(max(abs(term)) > tol && (n <- n+1) < maxit) {
> 	      term <- term*A*B/C * z/n
> 	      partial.sum <- partial.sum + term
> 	      A <- A+1
> 	      B <- B+1
> 	      C <- C+1
> 	    }
> 	    if(n == maxit) warning("convergence may not have been achieved")
> 	    partial.sum
> 	} 
> 


thank you very much for this: I guess the answer is to build
vectorization in, rather than to bolt it on afterwards!

The only advantage of my original function is that in mine, each
evaluation stops as soon as possible, whereas in Bill's improved
vectorized version, evaluation of each element continues until the
worst-case element finishes.  I imagine that this isn't much of a
problem, unless one deliberatly finds some arguments to make it so.

(incidentally, I am particularly taken with "term[] <- 1",
considerably better than ,y term <- 1+term*0 )-:

Perhaps I should have rephrased my question.  For a function which
absolutely cannot be written in a vectorized manner, eg pparz() below,
what is the best way to emulate true vectorization a la hyper() above?

[qparz() is a quantile function and pparz() inverts it to give the
CDF; try pparz(1:20,c(10,0.1, -0.1)) ]



qparz <- function(p , params) {
  C <- params[1]
  l1 <- params[2]
  l2 <- params[3]
  if(p <=  0) {return(0)  }
  if(p >=  1) {return(Inf)}
C*p^l1*(1-p)^l2
}

pparz <- function(x , params) {
     f <- function(p,a){qparz(p=p,params=a[1:3])-a[4]}
     options(warn = -1)
     if(length(x)==1) {
       return(     uniroot(f,c(0,1),tol=1e-5,a=c(params,x))$root)
     } else {
       return(sapply(x,pparz,params))
     }
}





-- 

Robin Hankin, Lecturer,
School of Geography and Environmental Science
Tamaki Campus
Private Bag 92019 Auckland
New Zealand

r.hankin at auckland.ac.nz
tel 0064-9-373-7599 x6820; FAX 0064-9-373-7042

as of: Thu Oct 24 08:37:00 NZDT 2002
This (linux) system up continuously for:  420 days, 14 hours, 19 minutes
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._