Integration in R

Please reply on list.

Hi David,
x[2] is the second variable, x2. It comes from the condition  
0<x1<x2<7.
No, it doesn't come from those conditions. It is being grabbed from  
some "x"-named object that exists in your workspace.

If your limits were 7 in both dimensions, then the code should be:

adaptIntegrate(f, lowerLimit = c(0, 0), upperLimit = c(7,7))
#----
$integral
[1] 228.6667

(At this point I was trusting R's calculus abilities more than yours.  
I wasn't too trusting of mine either, and so tried seeing if Wolfram  
Alpha would accept this expression:
  integrate 2/3 (x+y) over 0< x<7,  0<y<7

; which it did and calculating the decimal expansion of the exact  
fraction:

 > 686/3
[1] 228.6667
 >
David.

>
> Thanks.
>
> On 8 January 2013 18:11, David Winsemius <dwinsemius at comcast.net>  
> wrote:
>
> On Jan 8, 2013, at 9:43 AM, Naser Jamil wrote:
>
> Hi R-users.
>
> I'm having difficulty with an integration in R via
> the package "cubature". I'm putting it with a simple example here.   
> I wish
> to integrate a function like:
> f(x1,x2)=2/3*(x1+x2) in the interval 0<x1<x2<7. To be sure I tried it
> by hand and got 114.33, but the following R code is giving me  
> 102.6667.
>
> -------------------------------------------------------------------
> library(cubature)
> f<-function(x) { 2/3 * (x[1] + x[2] ) }
> adaptIntegrate(f, lowerLimit = c(0, 0), upperLimit = c(x[2],7))
>
>
> What is x[2]?  On my machine it was 0.0761, so I obviously got a  
> different answer.
>
> -- 
> David Winsemius, MD
> Alameda, CA, USA
>
>

David Winsemius, MD
Alameda, CA, USA
An embedded and charset-unspecified text was scrubbed...
Name: not available
URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20130108/92379a88/attachment-0001.pl>

Thanks. But then how to implement condition like 0<x1<x2<7? I would  
be happy to know that.
Multiply the function by the conditional expression:

 > f<-function(x) { 2/3 * (x[1] + x[2] )*(x[1] < x[2]) }
 > adaptIntegrate(f, lowerLimit = c(0, 0), upperLimit = c(7,7))
$integral
[1] 114.3333
David.
>
> On 8 January 2013 18:41, David Winsemius <dwinsemius at comcast.net>  
> wrote:
> Please reply on list.
>
>
> On Jan 8, 2013, at 10:27 AM, Naser Jamil wrote:
>
> Hi David,
> x[2] is the second variable, x2. It comes from the condition  
> 0<x1<x2<7.
>
> No, it doesn't come from those conditions. It is being grabbed from  
> some "x"-named object that exists in your workspace.
>
> If your limits were 7 in both dimensions, then the code should be:
>
> adaptIntegrate(f, lowerLimit = c(0, 0), upperLimit = c(7,7))
> #----
> $integral
> [1] 228.6667
>
> (At this point I was trusting R's calculus abilities more than  
> yours. I wasn't too trusting of mine either, and so tried seeing if  
> Wolfram Alpha would accept this expression:
>  integrate 2/3 (x+y) over 0< x<7,  0<y<7
>
> ; which it did and calculating the decimal expansion of the exact  
> fraction:
>
> > 686/3
> [1] 228.6667
> >
>
> -- 
> David.
>
>
>
>
> Thanks.
>
> On 8 January 2013 18:11, David Winsemius <dwinsemius at comcast.net>  
> wrote:
>
> On Jan 8, 2013, at 9:43 AM, Naser Jamil wrote:
>
> Hi R-users.
>
> I'm having difficulty with an integration in R via
> the package "cubature". I'm putting it with a simple example here.   
> I wish
> to integrate a function like:
> f(x1,x2)=2/3*(x1+x2) in the interval 0<x1<x2<7. To be sure I tried it
> by hand and got 114.33, but the following R code is giving me  
> 102.6667.
>
> -------------------------------------------------------------------
> library(cubature)
> f<-function(x) { 2/3 * (x[1] + x[2] ) }
> adaptIntegrate(f, lowerLimit = c(0, 0), upperLimit = c(x[2],7))
>
>
> What is x[2]?  On my machine it was 0.0761, so I obviously got a  
> different answer.
>
> -- 
> David Winsemius, MD
> Alameda, CA, USA
>
>
>
> David Winsemius, MD
> Alameda, CA, USA
>
>

David Winsemius, MD
Alameda, CA, USA

Thanks. But then how to implement condition like 0<x1<x2<7? I would be
happy to know that.

David implemented the condition by multiplying by x[1]<x[2] which in a numeric context is 0 when x[1]<x[2] en 1 when x[1]>=x[2]. That is what your requirement does.
The condition 0<x1<x2<7 has been split into three separate parts

0<x1<7
0<x2<7
x1<x2

The first two can be converted to arguments of adaptIntegrate.
The last is inserted into the function definition effectively making the function return 0 when x1>x2 which is what your inequality implies.

Berend
On 8 January 2013 18:41, David Winsemius <dwinsemius at comcast.net> wrote:

Please reply on list.

On Jan 8, 2013, at 10:27 AM, Naser Jamil wrote:

Hi David,
x[2] is the second variable, x2. It comes from the condition 0<x1<x2<7.

No, it doesn't come from those conditions. It is being grabbed from some
"x"-named object that exists in your workspace.

If your limits were 7 in both dimensions, then the code should be:

adaptIntegrate(f, lowerLimit = c(0, 0), upperLimit = c(7,7))
#----
$integral
[1] 228.6667

(At this point I was trusting R's calculus abilities more than yours. I
wasn't too trusting of mine either, and so tried seeing if Wolfram Alpha
would accept this expression:
integrate 2/3 (x+y) over 0< x<7,  0<y<7

; which it did and calculating the decimal expansion of the exact fraction:

686/3
[1] 228.6667

--
David.

Thanks.

On 8 January 2013 18:11, David Winsemius <dwinsemius at comcast.net> wrote:

On Jan 8, 2013, at 9:43 AM, Naser Jamil wrote:

Hi R-users.

I'm having difficulty with an integration in R via
the package "cubature". I'm putting it with a simple example here.  I wish
to integrate a function like:
f(x1,x2)=2/3*(x1+x2) in the interval 0<x1<x2<7. To be sure I tried it
by hand and got 114.33, but the following R code is giving me 102.6667.

------------------------------**------------------------------**-------
library(cubature)
f<-function(x) { 2/3 * (x[1] + x[2] ) }
adaptIntegrate(f, lowerLimit = c(0, 0), upperLimit = c(x[2],7))

What is x[2]?  On my machine it was 0.0761, so I obviously got a
different answer.

--
David Winsemius, MD
Alameda, CA, USA

David Winsemius, MD
Alameda, CA, USA

	[[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 implemented the condition by multiplying by x[1]<x[2] which in a numeric context is 0 when x[1]<x[2] en 1 when x[1]>=x[2].
OOPS!! Reverse the 0 and the 1 in that sentence (TRUE becomes 1 and FALSE becomes 0)

Berend

On 08-01-2013, at 22:00, Berend Hasselman <bhh at xs4all.nl> wrote:

?...
David implemented the condition by multiplying by x[1]<x[2] which  
in a numeric context is 0 when x[1]<x[2] en 1 when x[1]>=x[2].
OOPS!! Reverse the 0 and the 1 in that sentence (TRUE becomes 1 and  
FALSE becomes 0)
I didn't worry about which direction was used for the inequality in  
this case because of consideration of symmetry. (Thinking about it  
now, I still think I managed to choose the correct direction.)   
Obviously such a cavalier attitude should not be carried over to the  
general situation.
David Winsemius, MD
Alameda, CA, USA

On Jan 8, 2013, at 1:07 PM, Berend Hasselman wrote:

On 08-01-2013, at 22:00, Berend Hasselman <bhh at xs4all.nl> wrote:

?...
David implemented the condition by multiplying by x[1]<x[2] which  
in a numeric context is 0 when x[1]<x[2] en 1 when x[1]>=x[2].
OOPS!! Reverse the 0 and the 1 in that sentence (TRUE becomes 1 and  
FALSE becomes 0)
I didn't worry about which direction was used for the inequality in  
this case because of consideration of symmetry. (Thinking about it  
now, I still think I managed to choose the correct direction.)   
Obviously such a cavalier attitude should not be carried over to the  
general situation.
I hit the send button by mistake. (I originally misread Berend's  
correction.) I'm only posting this to make clear that it was myself,  
and not Berend, that I was accusing of having been cavalier.
David Winsemius, MD
Alameda, CA, USA