From: Luke Tierney <luke@stat.umn.edu> Message-Id: <9804131554.AA12497@nokomis.stat.umn.edu> Subject: Re: R-beta: S Compatibility (again) To: p.dalgaard@biostat.ku.dk (Peter Dalgaard BSA) Date: Mon, 13 Apr 1998 10:54:14 -0500 (CDT) In-Reply-To: <x23efh6aaf.fsf@blueberry.kubism.ku.dk> from "Peter Dalgaard BSA" at Apr 13, 98 04:55:36 pm Peter Dalgaard BSA wrote: > > Peter Dalgaard BSA <p.dalgaard@biostat.ku.dk> writes: > > > > Here is a cute example of what can be done in S but not in R. > > > Make a function for the pdf of an order statistic. > > > > > > > pdf.order <- function(n, r, pfun, dfun) { > > > con <- round(exp(lgamma(n + 1) - lgamma(r) - lgamma(n - r + 1))) > > > substitute( > > > function(x) { > > > Fx <- p(x) > > > K*Fx^r1*(1 - Fx)^nr*f(x) > > > }, > > > list(p = substitute(pfun), f = substitute(dfun), > > > r1 = r-1, nr = n-r, K = con) > > > ) > > > } > > > > pdf.order(9, 5, pnorm, dnorm) > > > function(x) > > > { > > > Fx <- pnorm(x) > > > 630 * Fx^4 * (1 - Fx)^4 * dnorm(x) > > > } > > > > > > The substitute()s to get unevaluated arguments do work but the > > > one to modify the function definitely does not. substitute() is > > > a very different kind of function in R from what it is in S. > > > > I see the problem. I suspect that there's a workaround, though. Will > > look at it. > > [a couple of hours later...] > > > pdf.order > function (n, r, pfun, dfun) > { > con <- round(exp(lgamma(n + 1) - lgamma(r) - lgamma(n - > r + 1))) > eval(eval(expression(substitute(function(x) { > Fx <- p(x) > K * Fx^r1 * (1 - Fx)^nr * f(x) > })), list(p = substitute(pfun), f = substitute(dfun), > r1 = r - 1, nr = n - r, K = con)), .GlobalEnv) > } > > pdf.order(9, 5, pnorm, dnorm) > function (x) > { > Fx <- pnorm(x) > 630 * Fx^4 * (1 - Fx)^4 * dnorm(x) > } > > pdf.order(9, 5, pnorm, dnorm)(0) > [1] 0.981772 > > This could have been much easier if substitute allowed a named list > like eval() does. I see no reason why it shouldn't... > > The necessity of the outer eval is more fundamental: R is sometimes > more stringent in distinguishing between objects and expressions that > evaluate to them. So without this necessary eval (<hehe..>) you'd get > an expression (or more precisely, a "call" object) that evaluates to a > function, not the function itself: > > I.e. > > S: > > mode(substitute(function()a+b)) > [1] "function" > > R: > > mode(substitute(function()a+b)) > [1] "call" > > mode(eval(substitute(function()a+b))) > [1] "function" > > You also need to eval() it in .GlobalEnv {or sys.frame(sys.parent(2))} > , otherwise the result will carry with it the environment of > pdf.order, which is a direct consequence of R's scoping rules. > Maybe we should think about starting an obfuscated S/R code contest? :-) The need to make heavy use of substitute or eval is, in my view, not a strengh of a language. The right tool for this job, in my view at least, is lexical scoping and closures. R privides these directly, so a much simpler definition is pdf.order<- function (n, r, pfun, dfun) { con <- round(exp(lgamma(n + 1) - lgamma(r) - lgamma(n - r + 1))) function(x) { Fx <- pfun(x) con * Fx^(r - 1) * (1 - Fx)^(n - r) * dfun(x) } } > pdf.order(9, 5, pnorm, dnorm)(0) [1] 0.981772 The free variables con, n, r, pfun and dfun in the returned function refer to the variables in the defining environment (this is lexical scope). Since the free variables in the value function are not modified you can do this in S if you abstact out the closure creation operation into a function MC, pdf.order<- function (n, r, pfun, dfun) { con <- round(exp(lgamma(n + 1) - lgamma(r) - lgamma(n - r + 1))) MC(function(x) { Fx <- pfun(x) con * Fx^(r - 1) * (1 - Fx)^(n - r) * dfun(x) }, list("con"=con, "pfun"=pfun, "dfun"=dfun, "r"=r, "n"=n)) } > pdf.order(9, 5, pnorm, dnorm)(0) [1] 0.981772 This isn't as clean as in R but it is clearer in its intent, to me at least, than the eval/substitute stuff. The definition of MC (make closure) I use is not based on substitute. Substitute, no matter how sophisticated, is in my view fundamentally the wrong tool for this job since it does not understand syntax. (It is a reasonable tool in this case since the function is very small, but in general it is not). Both S and R have syntactic constructs that can result in shadowing of variables -- substitute cannot understand these and will mess them up. An alternative is to add bindings for free variables as additional arguments to the function; in this example this produces > pdf.order(9, 5, pnorm, dnorm) function(x, con = 630, pfun = function(q, mean = 0, sd = 1) { ... }, dfun = function(x, mean = 0, sd = 1) { ... }, r = 5, n = 9) { Fx <- pfun(x) con * Fx^(r - 1) * (1 - Fx)^(n - r) * dfun(x) } This isn't perfect (real lexical scope is much better) but, in my opinion, it is a better approach than using substitute. Here is the definition of MC I use to implement this. It is fairly convoluted, but once you have it and understand conceptually what it does, which is quite simple, you don't need to look at it again. > MC function(f, env = NULL) { env <- as.list(env) if(mode(f) != "function") stop(paste("not a function:", f)) if(length(env) > 0 && any(names(env) == "")) stop(paste("all arguments are not named:", env)) fargs <- if(length(f) > 1) f[1:(length(f) - 1)] else NULL fbody <- f[length(f)] cf <- c(fargs, env, fbody) mode(cf) <- "function" return(cf) } Just my 2c worth :-) luke -- Luke Tierney University of Minnesota Phone: 612-625-7843 School of Statistics Fax: 612-624-8868 206 Church Street email: luke@stat.umn.edu Minneapolis, MN 55455 USA WWW: http://www.stat.umn.edu -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- 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@stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._