Question

I would like to use R to verify complicated derivative calculations. Here is my best attempt so far on a simple function:

f <- expression(a*log(x^2))
df.dx <- deriv(f, 'x')
(df.dx)
df.dx.byHand <- expression(2*a/x) # The derivative of f calculated by hand
(df.dx.byHand)
all.equal(df.dx, df.dx.byHand) 

The output of the above is

> (df.dx)
 expression({
.expr1 <- x^2
.value <- a * log(.expr1)
.grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
.grad[, "x"] <- a * (2 * x/.expr1)
attr(.value, "gradient") <- .grad
.value

> (df.dx.byHand)
expression(2 * a/x)

[1] "Component 1: target, current do not match when deparsed"

I have already checked

Explicit formula versus symbolic derivatives in R,

http://stat.ethz.ch/R-manual/R-patched/library/stats/html/deriv.html

Symbolic derivatives and simplification in R

and also

http://stat.ethz.ch/R-manual/R-patched/library/base/html/expression.html

since the core issue seems to be that I don't know how to extract from 'df.dx' the 'grad[, 'x']' part of the expression.

Due to IT issues, neither Sage nor YACAS are available to me.

Many thanks for suggestions!

Était-ce utile?

La solution

The problem in the question is not just extracting the expression but also the fact that the expression needs to be simplified.

1) If the "IT issues" mean you had trouble installing Ryacas see the Ryacas home page for troubleshotting tips and try this:

library(Ryacas)
x <- Sym("x")
a <- Sym("a")
identical(Eval(Simplify(deriv(a*log(x^2), x))), 2 * a / x) # TRUE

If you get FALSE then they still might be the same if the simplification does not reduce it to the same form as the hand produced one so you might want to try the next method in any case.

2) If the "IT issues" mean you are not allowed to install packages then just compare them on a grid:

df.dx.byHand <- function(x) 2*a/x

df.dx <- function(x) {}
body(df.dx) <- D(f,'x')

a <- 3
all.equal(df.dx.byHand(1:100), df.dx(1:100)) # TRUE

UPDATED Provided Ryacas solution

Licencié sous: CC-BY-SA avec attribution
Non affilié à StackOverflow
scroll top