IRTFM
IRTFM

Reputation: 263421

Integrate a function that has a function as a parameter in R

So I just answered a question and as soon as it was answered (correctly I think), the questioner deleted it. So here it is again:


I am new to R and need help getting this function to work. I need to create a function that can find the log of an mgf for any function and return values for the specified t. I've done lots of research and I found a lot of stuff telling me to use Vectorize() and to make sure I'm defining my parameters properly but i still can't seem to get this to work. I would love if someone could help me out!

I need to write a function that returns a numeric vector for natural log of mgf

# I'm using the expression 2*x as an example
# You can use any integrand as long as it is a function of x 
logmgf <- function(integrand, upper, lower, t) {
  expression <- function(x, t) {
   integrand * exp(x * t)
  }
  integrate <- integrate(expression(x, t), upper, lower)
  logmgf <- log(Vectorize(integrate[1]))
  return(logmgf)
}
logmgf(2 * x, upper = Inf, lower = 0, t = 0)

asked 2 hours ago xiong lui

Upvotes: 2

Views: 135

Answers (1)

IRTFM
IRTFM

Reputation: 263421

Let's try something more statistically or mathematically sensible, such as an un-normalized Normal distribution, namely the expression: exp(-x^2)

You are trying to create a new expression (actually an R "call") which will be parsed as the product of that expression times exp(x*t), so you need to a) deliver the function a real R language object and b) work with it using functions which will not mangle it. The quote-function will construct an expression that substitute can manipulate at the "language level". The function-function unfortunately will evaluate the "body" argument in a manner that will not honor your symbolic intent, so you need to use the body<- (function) which expects an expression on the right-hand side of the assignment operator. I'm going to leave in some of the debugging print(.) calls that I used to understand where I was going wrong in my earlier efforts:

logmgf <- function(integrand, upper, lower, t) {  
      expr <- substitute( integrand *exp(x*t), list(integrand=integrand) )
print(expr)  
      func <-  function(x ){}  # builds an empty function in x
      body(func)<- expr     # could have also set an environment
      # but in this case using envir=parent.frame() is not appropriate
 print(func)
      integral <- integrate( func, upper=upper, 
            # notice need to name the parameters
                            lower=lower 
            # else they would be positionally matched
            # (and therefore reversed in this case)
                           )$value 
            # the integrate fn returns a loist and the numeric part is in $value
      logmgf <- log(integral)
}
res <- logmgf(quote(exp(-x^2)), upper = Inf, lower = -Inf, t = 0)

> res
[1] 0.5723649

MGF's are integrated from -Inf to Inf (or for functions with restricted domains only over the x's with defined values).

I wanted to check that I would get the correct answer from a known argument so I added back the proper normalizing constant for a Normal distribution:

 mgf <- function(integrand, upper, lower, t) {  
       expr <- substitute( integrand *exp(x*t), list(integrand=integrand) )
       func <-  function(x ){}; body(func)<- expr
       integral <- integrate( func, upper=upper, lower=lower)$value
 }
 res <- mgf(quote((1/sqrt(2*pi))*exp(-x^2/2)), upper = Inf, lower = -Inf, t = 0)
 res
#[1] 1

Upvotes: 2

Related Questions