[R] Using uniroot() with output from deriv() or D()

[Gabor Grothendieck]

On Sat, Feb 26, 2011 at 2:16 AM, jjheath <heath_jeremy at hotmail.com> wrote:
> I need to find the root of the second derivative of many curves and do not
> want
> to cut and paste the expression results from the deriv() or D() functions
> every time.  Below is an
> example.  What I need to do is refer to "fn2nd" in the uniroot() function,
> but when I
> try something like uniroot(fn2nd,c(0,1)) I get an error, so I have resorted
> to pasting
> in the expression, but this is highly inefficient.
>
> Thanks,  J
>
> y <- c(9.9,10,10,9.5,9,6,3,1,0,0,0)
> b <- seq(0,1,by=0.1)
> dat <- data.frame(y = y, b = b)
> plot(y ~ b, xlim = c(0,1), ylim= c(-12,12))
> fn <- nls(y ~ asym/(1 + exp(p * b - q)),data = dat, list(asym=10,p=16,q=7),
> trace=T)
> curve(10.001/(1 + exp(14.094 * x - 7.551)), from = 0, to = 1, add = T)
> fn2 <- expression(10.001/(1 + exp(14.094 * x - 7.551)))
> fn2nd <- D(D(fn2, "x"), "x")
> ex <- seq.int(from=0, to=1, length.out = 1000)
> y1 <- eval(fn2nd, envir = list(x = ex), enclos = parent.frame())
> lines(ex, y1, type = "l")
> r <- uniroot(function(x) -(10.001 * (exp(14.094 * x - 7.551) * 14.094 *
> 14.094)/(1 + exp(14.094 *
>    x - 7.551))^2 - 10.001 * (exp(14.094 * x - 7.551) * 14.094) *
>    (2 * (exp(14.094 * x - 7.551) * 14.094 * (1 + exp(14.094 *
>        x - 7.551))))/((1 + exp(14.094 * x -
> 7.551))^2)^2),interval=c(0,1),tol = 0.0001)
> r$root
> abline(h=0, col = "red")
> abline(v=r$root, col = "green")
> arrows(0.6, -2, r$root,0, length = 0.1, angle = 30, code = 2, col = "red")
> text(0.765,-2.3,paste("b = ",r$root,sep=""))

Try this:

f <- function(x) {}
body(f) <- fn2nd

-- 
Statistics & Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.com