[R] Gentleman and Ihaka , 2000 paper question

[Robert Gentleman]

not sure just what you want, but here are some snippets

newton <-
   function(lfun, est, tol = 1e-7, niter = 500) {
     cscore <- lfun$score(est)
     if (abs(cscore) < tol)
       return(est)
     for (i in 1:niter) {
       new <- est - cscore / lfun$d2(est)
       cscore <- lfun$score(new)
       if (abs(cscore) < tol)
         return(new)
       est <- new
     }
     stop("exceeded allowed number of iterations")
   }

with the likelihood function

Rmklike <-
   function(data) {
     n <- length(data)
     sumx <- sum(data)
     lfun <- function(mu) n * log(mu) - mu * sumx
     score <- function(mu) n / mu - sumx
     d2 <- function(mu) -n / mu^2
     list(lfun = lfun, score = score, d2 = d2)
   }


Leeds, Mark (IED) wrote:
> In their paper, "Lexical Scope and Statistical Computing", the authors (
> Gentleman and Ihaka ) go to great length explaining why R's use of
> lexical scoping creates advantages when doing statistical computations.
> 
> If anyone has or is familiar with this paper, could they provide the
> main program code for how the "newton" function would be called in their
> example on page 500 of the paper.  The authors are extremely clear in
> their writing and the paper is quite an eye opener for me but it seems
> like lfun somehow needs to be initialized so that it "grabs" the
> environment of "Rmklike". 

   Rmlike is called to create the likelihood function, and since that 
function is defined in the body of Rmlike, it has the evaluation 
environment by default. And any return value will have the correct 
environment.

> I'm not sure how one would go about doing this so I am wondering what
> the main program that calls "newton" would be
> if there was one. Thanks.
> 

So,

   data=rexp(10, rate=.3)
   lf = Rmklike(data)
   newton(lf, .1)


seems to do the trick (and surprisingly still works, as written some 
unbelievably long time ago)

  Robert

> 
> 	
> Mark
> --------------------------------------------------------
> 
> This is not an offer (or solicitation of an offer) to buy/se...{{dropped}}
> 
> ______________________________________________
> R-help at stat.math.ethz.ch 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.
> 

-- 
Robert Gentleman, PhD
Program in Computational Biology
Division of Public Health Sciences
Fred Hutchinson Cancer Research Center
1100 Fairview Ave. N, M2-B876
PO Box 19024
Seattle, Washington 98109-1024
206-667-7700
rgentlem at fhcrc.org