[R] Questions about lexical scope
Lu Chi-Hsien Joseph
cjlu at ibm.stat.ncku.edu.tw
Tue Aug 6 17:28:28 CEST 2002
Dear R-users,
The numerical integration example given in Gentleman and Ihaka (2000),
"Lexical Scope and Statistical Computing," JCGS, 9, 491-508,
is very interesting and helpful in understanding how lexical scope
is about.
However, I got some questions that I just can't figure out.
First all, allow me to copy the two functions given by the authors:
midpoint <- function(f, a, b, n=100) {
h <- (b - a)/n
(b - a)*mean(sapply(seq(a + h/2, b - h/2, length=n), f))
}
bv.integrate <- function(f, a, b, n=100, rule=midpoint) {
g <- function(y) {
fx <- function(x) f(x, y)
rule(fx, a[2], b[2], n)
}
rule(g, a[1], b[1])
}
I modified the function name from "integrate" to "bv.integrate".
To integrate f(x,y)=sin(x + y) over the unit square,
it works:
> bv.integrate(function(x, y) sin(x + y), c(0,0), c(1,1))
[1] 0.773651
(compare with .7736445 calculated from 2sin(1)*(1 - cos(1)))
Then I write a function using trapezoidal rule as follows:
trapezoidal <- function(f, a, b, n=100) {
h <- (b - a)/n
i <- 1:n
sum(f(a + (i - 1)*h), f(a + i*h))*h/2
}
I checked with
> trapezoidal(sin, 0, 1)
[1] 0.4596939
> 1 - cos(1)
[1] 0.4596977
and
> trapezoidal(dnorm, -3, 3)
[1] 0.9972922
> diff(pnorm(c(-3, 3)))
[1] 0.9973002
Then, this is what I got by plugged in "trapezoidal" for the "rule"
in bv.integrate():
> bv.integrate(function(x, y) sin(x + y), c(0,0), c(1,1), rule=trapezoidal)
[1] 0.007080675
Hundred time smaller!
Answer can be "improved" by increasing n, but always 100 times smaller.
Anything wrong with the way I wrote trapezoidal()?
I just can't figure out any difference, in terms of input/output,
between midpoint() and trapezoidal().
Then I tried to use integrate() as the rule in bv.integrate(),
by defining a wrap-up function:
intg <- function(f, a, b, n=100) integrate(f, a, b, n)$value
to have the value as the only output.
But, I got
> bv.integrate(function(x, y) sin(x + y), c(0,0), c(1,1), rule=intg)
Error in integrate(f, a, b, n) : evaluation of function gave a result
of wrong length
I tried to use debugger() to find out what's wrong,
but not succeeded.
There must be something I missed in understanding the lexical scope,
however, what's going on with the function I wrote?
What should be the way to correct them?
Certainly, I know integrate() should be the one to use.
These are questions from my preparation of programming exercise
of using R in my statistical computation class.
I'll be very appreciated for anyone who might help me on these
questions.
Best regards,
C. Joseph Lu
Department of Statistics
National Cheng-Kung University
Tainan, Taiwan, ROC
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
More information about the R-help
mailing list