Reputation: 263421
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
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