CODEWITHSUNDEEP
rderivativequotation-marks
PlaymoBill
PlaymoBill

Reputation: 23

Computing the nth derivative of a function

To compute a probability, I have to compute derivatives (and then evaluate) like $\frac{\partial^5 f}{\partial x_1^2 \partial x_2^3}$ where $f$ is a polynomial function. The problem is that the order of the derivative is likely to vary as well as the list of variables with respect to which the derivative is computed.

I already tried with rSymPy and Ryacas and it works... until the number of variables becomes to important. So I have to look for a different solution. I tried with the DD() function indicated in the documentation of deriv() and using this function iteratively seems to be fine (and unexpectedly more efficient than with rSymPy and Ryacas).

My problem is to create the DD(DD(DD(...my.expr...,"xi",ni),"xj",nj),"xk",nk) command. I tried the following code:

step1 <- function(k) paste0(",x", k, ",", r[k]-1, ")", collapse="")
step2 <- function(expr) {
           paste0(paste0(rep.int("DD(",u), collapse=""), expr, 
           paste0(sapply(t,f4), collapse=""), collapse="") }
step2(f)

where r is a vector indicating the order of derivation for each variable, t a subset of that vector, u <- length(t) and f is an expression object. This solution does not works because quotation marks are missing around variable names. Indeed I get for instance (I dropped the function from the code):

DD(DD(DD(DD(DD(my.expr,x1,1),x7,1),x9,2),x10,1),x11,1)

instead of:

DD(DD(DD(DD(DD(my.expr,"x1",1),"x7",1),"x9",2),"x10",1),"x11",1)

I tried adding \" in my function step1, but I have then a problem with the computation of the derivative. Any suggestion to fix this problem?

PS: it would surely be easier with a loop, but I would like to avoid if possible.
PS2: Sorry for LaTeX code.

Upvotes: 2

Views: 841

Answers (1)

Ben Bolker
Ben Bolker

Reputation: 226732

I think this extension works. The trick is to not start going back and forth between expressions and strings ...

DD <- function(expr, names, order = 1, debug=FALSE) {
    if (any(order>=1)) {  ## do we need to do any more work?
        w <- which(order>=1)[1]  ## find a derivative to compute
        if (debug) {
            cat(names,order,w,"\n")
        }
        ## update order
        order[w] <- order[w]-1
        ## recurse ...
        return(DD(D(expr,names[w]), names, order, debug))
    }
    return(expr)
}

Some tests:

DD(expression(x^2*y^3+z),c("x","y"),c(1,1))
## 2 * x * (3 * y^2)
DD(expression(x^2*y^3+z),c("x","y"),c(2,1))
## 2*3*(y^2)
DD(expression(x^2*y^3+z),c("x","y"),c(2,2))
## 2*(3*(2*y))
DD(expression(x^2*y^3+z),c("x","y"),c(2,3))
## 2*(3*2)
DD(expression(x^2*y^3+z),c("x","y"),c(2,4))
## 0

I hadn't noticed previously that you were differentiating a polynomial -- in that special case there's a much simpler answer (hint, represent the polynomial as a sequence of vectors that give the coefficients of orders of different terms). But you may not need that efficient an answer ...

Upvotes: 2

Related Questions