2

In Octave, I use the following to differentiate a function of 2-variables and then substitute 0 for x using subs(). Basically in doing moment-generating-function, Taylor series expansion, etc, we differentiate and then substitute some a for x. I am not able to find the equivalent substitution function in R. Can you please let me know how to do it? Thanks

pkg load symbolic; #octave symbolic package
syms lamb, x; #declare lamb, x symbols
mgf = lamb / (lamb - x); #moment generating function of exponential
mgf1 = diff(mgf, x, 1); #1st differivative
mgf1_0 = subs(mgf1, x, 0); #substitute 0 for E(X)
Martijn Pieters
  • 1,048,767
  • 296
  • 4,058
  • 3,343
limestreetlab
  • 173
  • 1
  • 11
  • your example would be a little more useful for future readers if you include the necessary package-loading and variable-defining steps ... – Ben Bolker Jul 04 '20 at 23:47
  • added pkg. look fwd to your future helps on R, which I am rather confused about in many aspects. – limestreetlab Jul 05 '20 at 08:08

3 Answers3

4

Using base R:

f <- quote( lambda / (lambda - x) )
Df <-  D(f, "x")

do.call("substitute", list(Df, list(x = 0)))
## lambda/(lambda - 0)^2

or we can evaluate Df at specific x and lambda values:

eval(Df, list(x = 0, lambda = 3))
## [1] 0.3333333

Create function to represent result

We can use Df to define an R function der which evaluates the derivative at specific x and lambda values.

der <- function(x, lambda) {}
body(der) <- Df
der(0, 3)
## [1] 0.3333333

Currying

If we want to set x to 0 and create the resulting function of lambda

make_der0 <- function(x = 0) function(lambda) der(x, lambda)
der0 <- make_der0()
der0(3)
## [1] 0.3333333

This is known as currying or partialling and various packages have this as well:

library(functional)
der0a <- Curry(der, x = 0)
der0a(3)
## [1] 0.3333333

library(purrr)
der0b <- partial(der, x = 0)
der0b(3)
## [1] 0.3333333
G. Grothendieck
  • 254,981
  • 17
  • 203
  • 341
  • Cool answer to use `quote` and `substitute`, +1! – ThomasIsCoding Jul 03 '20 at 18:12
  • This is nice but doesn't simplify the result. (I find it a little bit ironic that you wrote `rSymPy`, which my answer uses ...) – Ben Bolker Jul 03 '20 at 18:17
  • by the way, I'm a little confused about the development status of `rSymPy`. I posted an issue at https://github.com/cjgb/rSymPy but I'm not sure if it's current: http://rsympy.googlecode.com looks broken ... ? – Ben Bolker Jul 03 '20 at 18:59
  • I'd be willing, but I'd only be able to do dumb R-level maintenance stuff (e.g., maintaining CRAN compliance). If anything broke in the internals the probability I'd be able to fix it would be <0.001 ... – Ben Bolker Jul 03 '20 at 20:07
  • @BB, Have moved discussion to email. – G. Grothendieck Jul 03 '20 at 22:20
  • I use eval(f,values) when substituting specific values for all symbols, but many function-to-function problems require substituting for partial symbols only. As an everyday user, R seems to make me think about how to do something as common as symbol substitution, while I can just think about what to do in Octave. – limestreetlab Jul 05 '20 at 08:13
2

This appears to be using the symbolic package, which in turn depends on the sympy Python library. R doesn't have built-in symbolic manipulation capabilities, but it does (TIL) have an rSymPy package that works similarly.

## https://kevinkotze.github.io/mm-tut1-symbolic/
library(rSymPy)

rSymPy doesn't have a built-in subs() function, so we'll define one:

subs <- function(expr,x,y) {
    Sym(expr,".subs(",x,",",y,")")
}

Also useful to define this (not sure if there's a better way):

r_eval <- function(s,eval_list) {
   eval(parse(text=sympy(unclass(s))), eval_list)
}

The rest of the code looks almost identical to the Octave code. Define variables:

## NOTE: 'lambda' is a reserved word in Python, so trying to use it as 
## a variable gives rise to confusing errors ...
lam <- Var("lam")
x <- Var("x")
f <- lam/(lam-x)

Differentiate:

mgf1 <- deriv(f,x,1)

Substitute:

subs(mgf1,x,0)

There are other interfaces from R to symbolic math engines, e.g. Ryacas.

If you want to compute the second moment by calculating f''(0), that's not much harder:

v <- subs(deriv(f,x,2),x,0)  ## "2/lam**2"
r_eval(v, list(lam=5))  ## 0.08
Ben Bolker
  • 211,554
  • 25
  • 370
  • 453
  • I have been playing with moment generating functions and taylor series expansions etc in Octave, and always relied on subs(f,x,value). I didn't expect it is such a hassle in R for a seemingly simple operation. – limestreetlab Jul 03 '20 at 18:17
  • 1
    Well, neither platform has native symbolic computation, and I've just showed how you can do it just about as easily in R as in Octave ... I should have defined the `subs()` function first, then the code would have been pretty much identical to yours (except that you left out the package-loading and variable-definition steps ... – Ben Bolker Jul 03 '20 at 18:20
  • A similar issue when doing college math in R or Octave, can you please help shed some light on https://stackoverflow.com/questions/62853882/how-to-calculate-a-bayes-estimator-using-octave-or-matlab – limestreetlab Jul 12 '20 at 02:19
0

Maybe you can try package Deriv

library(Deriv)
f <- function(l,x) l/(l-x)
fprim <- Deriv(f,"x")(2,0)

such that

> fprim(2,0)
[1] 0.5
ThomasIsCoding
  • 96,636
  • 9
  • 24
  • 81
  • didn't know that package. Can it substitute x only? I see you substituting l=2, x=0, which can be done using eval(mgf1). – limestreetlab Jul 04 '20 at 01:15