[R] Taking Integral and Optimization using Integrate, Optim and maxNR
Berend Hasselman
bhh at xs4all.nl
Tue Jun 7 20:40:34 CEST 2011
Dear All, Hello!
I have some questoins in R programming as follows:
Question 1- How to take the integral of this function with respect to y,
such that x would appear in the output after taking integral.
f(x,y)=(0.1766*exp(-exp(y+lnx))*-exp(y+lnx))/(1-exp(-exp(y+lnx))) y in
(-6.907,-1.246)
It is doable in maple but not in R. At least I could not find the way.
p.s: result from maple is:
g(x)=dilog*exp(0.001000755564*x)+0.5*ln(exp(0.001000755564*x))^2-dilog*(exp(0.2876531111*x))-0.5*ln(exp(0.2876531111*x))^2
Where dilog=integral(log(t)/(1-t)) for t in (1,x)
Question 2- Then I want to optimize (maximize) the answer of the above
integral for x. Assum x in (0,100)
The answer should be something between 26 and 27.
What I did is: got answer of integral(f(x,y)) from maple which is g(x), and
then applied it in R. The code is as follows:
##In the case n=1
library(stats)
require(graphics)
integrand=function(t){(log(t))/(1-t)}
#dilog=integrate(integrand, lower=1, upper=x)
fr <- function(x){(integrate(integrand, lower=1,
upper=x)$value)*(exp(0.001000755564*x))+0.5*log(exp(0.001000755564*x))^2-(integrate(integrand,
lower=1,
upper=x)$value)*(exp(0.2876531111*x))-0.5*log(exp(0.2876531111*x))^2}
optim(20, fr, NULL, method = "BFGS")
Question 2-1-Default by optim is to minimize a function, how I can use it to
maximization?
Question 2-2- The above code faced with errors, and did not work. What I
guess is there is something wrong with taking integral. The output of
integrate function in R is some sort of thing. I had to somehow tell it,
just take the value and forget about the others. But I guess it is still
something wrong with it. I also tried maxNR, but is didn't work either.
/quote>
The integral of dilog (should) exist (according to the books).
dilog <- function(t) log(t)/(1-t)
integrate(dilog,1,2)
-0.822467 with absolute error < 9.1e-15
Now if you do:
z <- integrate(dilog,1.000,2)
str(z)
you will see that z is a list:
List of 5
$ value : num -0.822
$ abs.error : num 9.13e-15
$ subdivisions: int 1
$ message : chr "OK"
$ call : language integrate(f = dilog, lower = 1, upper = 2)
- attr(*, "class")= chr "integrate"
If you need the result of integrate then it is standard R to use z$value to
get the value of the integral.
Berend
--
View this message in context: http://r.789695.n4.nabble.com/Taking-Integral-and-Optimization-using-Integrate-Optim-and-maxNR-tp3577923p3580456.html
Sent from the R help mailing list archive at Nabble.com.
More information about the R-help
mailing list