An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20131224/55e4c784/attachment.pl>
Simplifying an expression with an integral
2 messages · Aurélien Philippot, David Winsemius
On Dec 24, 2013, at 9:39 AM, Aur?lien Philippot wrote:
Dear R experts,
I am computing the following integral.
[image: \int_{1,100} \frac{1}{x+Max(x-50,0)} g(x)dx], where g is the
density of the standard normal, and [1,100] is the domain.
1) I use the following code which works fine:
integrand1<- function(x){1/x*dnorm(x)}
integrand2<- function(x){1/(2*x-50)*dnorm(x)}
res1<- integrate(integrand1, lower=1, upper=50)$value+
integrate(integrand2, lower=50, upper=100)$value
res1
0.1116
In other words, I split the max function depending on the value of x in the
domain.
2) Alternatively, I can also compute it by vectorizing the max function
integrand<- function (x) ifelse(x<50, 1/x*dnorm(x) , 1/(2*x-50)*dnorm(x))
res4<- integrate(integrand, lower=1, upper=100)$value
res4
0.1116
3) However, in both cases, the syntax is a little bit heavy and not very
convenient if I want to add more integrals, all of them with a max in the
integrand. Is there a way to have a more concise syntax?
Using max or pmax directly in the definition of the integrand does not work.
For example:
integrand<- function(x) {1/(x+max(x-50,0)*dnorm(x))}
Is boolean math more to your liking?
integr_both<- function(x){
(x < 50)*(1/x*dnorm(x))+
(x>= 50 & x<100)*( 1/(2*x-50)*dnorm(x) ) }
res_b<- integrate(integr_both, lower=1, upper=100)$value
res_b
[1] 0.1116587
David > > res<- integrate(integrand, lower=1, upper=100)$value > res > 4.60517 > > Thank you for any suggestion, and merry Christmas! > > [[alternative HTML version deleted]] > > ______________________________________________ > R-help at r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. David Winsemius Alameda, CA, USA