[R] getting variable length numerical gradient

Antonio, Fabio Di Narzo antonio.fabio at gmail.com
Tue Sep 27 15:27:38 CEST 2005


2005/9/26, Dimitris Rizopoulos <dimitris.rizopoulos at med.kuleuven.be>:
> AFAIK the deriv() is for symbolic derivatives, so I don't know if it
> will work in your case.

numericDeriv surely computes numerical gradient. The problem here is
its interface, which requires the vector of symbols involved in the
expression to differenziate.


>
> Best,
> Dimitris
>
>
> ----- Original Message -----
> From: "Antonio, Fabio Di Narzo" <antonio.fabio at gmail.com>
> To: "Dimitris Rizopoulos" <dimitris.rizopoulos at med.kuleuven.be>
> Cc: <R-help at stat.math.ethz.ch>
> Sent: Monday, September 26, 2005 10:58 AM
> Subject: Re: getting variable length numerical gradient
>
>
> Tnx very much Dimitris,
> your code does what I need. I've just adapted it to my needs (e.g., I
> don't deal with scalar functions), and so solved my problem.
>
> Given this, is there a way to use the deriv function in the base
> package, within this context (variable length vector of indipendent
> variables)?
>
> Best,
> Antonio, Fabio Di Narzo.
>
> On 9/25/05, Dimitris Rizopoulos <dimitris.rizopoulos at med.kuleuven.be>
> wrote:
> > maybe you can find the following function useful (any comments are
> > greatly appreciated):
> >
> > fd <- function(x, f, scalar = TRUE, ..., eps =
> > sqrt(.Machine$double.neg.eps)){
> >     f <- match.fun(f)
> >     out <- if(scalar){
> >         if(length(f0 <- f(x, ...)) != length(x))
> >             stop("'f' must be vectorized")
> >         x. <- x + eps * pmax(abs(x), 1)
> >         c(f(x., ...) - f0) / (x. - x)
> >     } else{
> >         n <- length(x)
> >         res <- array(0, c(n, n))
> >         f0 <- f(x, ...)
> >         ex <- pmax(abs(x), 1)
> >         for(i in 1:n){
> >             x. <- x
> >             x.[i] <- x[i] + eps * ex[i]
> >             res[, i] <- c(f(x., ...) - f0) / (x.[i] - x[i])
> >         }
> >         res
> >     }
> >     out
> > }
> >
> >
> > ## Examples
> >
> > x <- seq(-3.3, 3.3, 0.1)
> > all.equal(fd(x, pnorm, mean = 0.5), dnorm(x, mean = 0.5))
> >
> >
> > # Approximate the Hessian matrix for a logistic regression
> >
> > # the score vector function
> > gn <- function(b, y, X){
> >     p <- as.vector(plogis(X %*% b))
> >     -colSums(X * (y - p))
> > }
> >
> > # We simulate some data and fit the logistic regression
> > n <- 800
> > x1 <- runif(n,-3, 3); x2 <- runif(n, -3, 3)
> > pr <- plogis(0.8 + 0.4 * x1 - 0.3 * x2)
> > y <- rbinom(n, 1, pr)
> > fm <- glm(y ~ x1 + x2, binomial)
> >
> > ## The Hessian using forward difference approximation
> > fd(fm$coef, gn, scalar = FALSE, y = y, X = cbind(1, x1, x2))
> >
> > ## The true Hessian
> > solve(summary(fm)$cov.unscaled)
> >
> >
> > I hope it helps.
> >
> > Best,
> > Dimitris
> >
> > ----
> > Dimitris Rizopoulos
> > Ph.D. Student
> > Biostatistical Centre
> > School of Public Health
> > Catholic University of Leuven
> >
> > Address: Kapucijnenvoer 35, Leuven, Belgium
> > Tel: +32/(0)16/336899
> > Fax: +32/(0)16/337015
> > Web: http://www.med.kuleuven.be/biostat/
> >      http://www.student.kuleuven.be/~m0390867/dimitris.htm
> >
> >
> > ----- Original Message -----
> > From: "Antonio, Fabio Di Narzo" <antonio.fabio at gmail.com>
> > To: <R-help at stat.math.ethz.ch>
> > Sent: Sunday, September 25, 2005 11:37 AM
> > Subject: [R] getting variable length numerical gradient
> >
> >
> > > Hi all.
> > > I have a numerical function f(x), with x being a vector of generic
> > > size (say k=4), and I wanna take the numerically computed
> > > gradient,
> > > using deriv or numericDeriv (or something else).
> > >
> > > My difficulties here are that in deriv and numericDeric the
> > > function
> > > is passed as an expression, and one have to pass the list of
> > > variables
> > > involved as a char vector... So, it's a pure R programming
> > > question.
> > >
> > >
> > > Have a nice sunday,
> > > Antonio, Fabio Di Narzo.
> > >
> > > ______________________________________________
> > > 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
> > >
> >
> >
> > Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
> >
> >
>
>
> Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
>
>




More information about the R-help mailing list