RuleBasedIntegration / Rubi

Rubi for Mathematica
http://rulebasedintegration.org
MIT License
225 stars 22 forks source link

Integrals in the context of moments of a random variables (Rubi 4.16.0.4) #6

Closed hedgehog closed 5 years ago

hedgehog commented 5 years ago

Thanks for all your efforts into making Rubi so powerful.

I'm not sure if the following identifies an issue with Rubi or the keyboard operator ;)

I'm trying to verify some calculations with MMA, and I have been unable to get the same results as in the literature.
So I'm now trying with Rubi and have gone back to the source reference:

Heyde, CC. (1963), "On a property of the lognormal distribution", Journal of the Royal Statistical Society, Series B, 25 (2): 392–393 (the pdf here: http://link.springer.com/chapter/10.1007%2F978-1-4419-5823-5_6 )

The relevant excerpt is:

image

If I try this in Rubi....

f   =  1/((x - c) b Sqrt[2 \[Pi]]) Exp[-((Log[x - c] - a)^2/(2 b^2))]
Int[x*f, x]
Int[x^2*f, x]
Int[x^3*f, x]

.... I get results that are close. Acknowledge the Rubi does not accept restrictions on the domains of the parameters, nonetheless I'm not sure how to interpret the un-evaluated Int[...] terms:

-(1/2) c Erf[(a - Log[-c + x])/(Sqrt[2] b)] - 
 1/2 E^(a + b^2/2) Erf[(a + b^2 - Log[-c + x])/(Sqrt[2] b)]

-(1/2) c^2 Erf[(a - Log[-c + x])/(Sqrt[2] b)] - 
 1/2 c E^(a + b^2/2) Erf[(a + b^2 - Log[-c + x])/(Sqrt[2] b)] + 
 Int[E^(-((-a + Log[-c + x])^2/(2 b^2))) x, x]/(b Sqrt[2 \[Pi]])

-(1/2) c^3 Erf[(a - Log[-c + x])/(Sqrt[2] b)] - 
 1/2 c^2 E^(a + b^2/2) Erf[(a + b^2 - Log[-c + x])/(Sqrt[2] b)] + (
 c Int[E^(-((-a + Log[-c + x])^2/(2 b^2))) x, x])/(b Sqrt[2 \[Pi]]) + 
 Int[E^(-((-a + Log[-c + x])^2/(2 b^2))) x^2, x]/(b Sqrt[2 \[Pi]])

Is this the expected result?

I arrived at the above after finding my attempt to use the following did not return any result:

y = If[x <= c, 0, x]
Int[y^2 * f,x]
(*etc...*)

Then Rubi returns the result un-evaluated.

halirutan commented 5 years ago

The un-evaluated Int expression means that Rubi cannot integrate this part. We will check, why this is. In the mean-time, you can use the normal Integrate. When I tested it in 11.3, Mathematica was able to find antiderivatives:

In[22]:= f = 
  1/((x - c) b Sqrt[2 \[Pi]]) Exp[-((Log[x - c] - a)^2/(2 b^2))];
Integrate[x*f, x]
Integrate[x^2*f, x]
Integrate[x^3*f, x]

Out[23]= -(1/2) c Erf[(a - Log[-c + x])/(Sqrt[2] b)] - 
 1/2 E^(a + b^2/2) Erf[(a + b^2 - Log[-c + x])/(Sqrt[2] b)]

Out[24]= -(1/2) c^2 Erf[(a - Log[-c + x])/(Sqrt[2] b)] - 
 1/2 E^(a + b^2/
   2) (2 c Erf[(a + b^2 - Log[-c + x])/(Sqrt[2] b)] + 
    E^(a + (3 b^2)/2) Erf[(a + 2 b^2 - Log[-c + x])/(Sqrt[2] b)])

Out[25]= -(1/2) c^3 Erf[(a - Log[-c + x])/(Sqrt[2] b)] - 
 1/2 E^(-(a^2/(
   2 b^2))) (3 c^2 E^((a + b^2)^2/(2 b^2))
      Erf[(a + b^2 - Log[-c + x])/(Sqrt[2] b)] + 
    3 c E^((a + 2 b^2)^2/(2 b^2))
      Erf[(a + 2 b^2 - Log[-c + x])/(Sqrt[2] b)] + 
    E^((a + 3 b^2)^2/(2 b^2))
      Erf[(a + 3 b^2 - Log[-c + x])/(Sqrt[2] b)])
hedgehog commented 5 years ago

Thanks @halirutan. Just to clarify: MMA (mathStatica's Expect[...]) could not handle Heyde's density g(x). So I turned to Rubi to calculate it. I then reduced it to f(x) and found Rubi could not resolve/evaluate the reduced expression.

It may also be worth checking if Heyde's g(x) can be resolved once f(x) is able to be calculated?

To be clear I assume the un-evaluated Int[...] you refer to is from the Int[x*f,z] formulation.

Would you expect Int[If[x<=c,0,x]*f,x] to succeed?

Finally, I note the MMA 11.2 Integrate[x^6 * f, x] returns un-evaluated. Lower order powers seem to be evaluated.

halirutan commented 5 years ago

To be clear I assume the un-evaluated Int[...] you refer to is from the Int[x*f,z] formulation.

It does not really matter where it comes from. An unevaluated Int always indicates that Rubi cannot find a rule to compute the antiderivative for this expression. It might be that this is flaw in the rules and Rubi should be able to compute it but the bug prevents the finding of the correct rule or it might be that Rubi has no rule for this expression. I'll wait what Albert says.

Would you expect Int[If[x<=c,0,x]*f,x] to succeed?

No. As far as I can tell, Rubi has not knowledge of these constructs. Here, you are better of with Integrate if it can do it.

Finally, I note the MMA 11.2 Integrate[x^6 * f, x] returns un-evaluated.

MMA 11.3 seems to be able to solve this. Can you try a hybrid approach? First, you integrate with Rubi and then you feed the un-evaluated integrals to Mathematica's Integrate like this:

f = 1/((x - c) b Sqrt[2 \[Pi]]) Exp[-((Log[x - c] - a)^2/(2 b^2))];
Int[x^6*f, x]

%/. Defer[Int] -> Integrate
hedgehog commented 5 years ago

Wow...

%/. Defer[Int] -> Integrate

Great suggestion, I wasn't aware of that usage. Thanks!

AlbertRich commented 5 years ago

Thanks for your kind words about Rubi, and for reporting this deficiency in Rubi 4.16.0. I spent the day generalizing the exponential integration rules required to integrate expressions of which your examples are a special case. So the next version of Rubi will be able to integrate all expressions of the form

(g + h*x)^m * F^(f*(a + b*Log[c*(d + e*x)^n])^2)

if m is a nonnegative integer or e g-d h=0.

Then for your three integration problems Rubi will produce

(-(1/2))*c*Erf[(a - Log[-c + x])/(Sqrt[2]*b)] - (1/2)*E^(a + b^2/2)*
  Erf[(a + b^2 - Log[-c + x])/(Sqrt[2]*b)]

(-(1/2))*c^2*Erf[(a - Log[-c + x])/(Sqrt[2]*b)] - 
 c*E^(a + b^2/2)*Erf[(a + b^2 - Log[-c + x])/(Sqrt[2]*b)] - (1/2)*
  E^(2*(a + b^2))*Erf[(a + 2*b^2 - Log[-c + x])/(Sqrt[2]*b)]

(-(1/2))*c^3*Erf[(a - Log[-c + x])/(Sqrt[2]*b)] - (3/2)*c^2*
  E^(a + b^2/2)*Erf[(a + b^2 - Log[-c + x])/(Sqrt[2]*b)] - (3/2)*c*
  E^(2*(a + b^2))*Erf[(a + 2*b^2 - Log[-c + x])/(Sqrt[2]*b)] - (1/2)*
  E^(3*a + (9*b^2)/2)*Erf[(a + 3*b^2 - Log[-c + x])/(Sqrt[2]*b)]

Aloha from Hawaii, Albert

halirutan commented 5 years ago

This fix will be included in the next release.

hedgehog commented 5 years ago

Thanks @AlbertRich, much appreciated. You're a hard marker... 'deficiency'. I would have said logic gap ;)

What I am really curious about is whether, for g(x) defined in the snippet from Heyde's paper above, Rubi is able to calculate:

Integrate[x*g, x]
Integrate[x^2*g, x]
Integrate[x^3*g, x]

and whether they agree with the results from Integrate[(x^n)*f, x], n=1,2,3,...

AlbertRich commented 5 years ago

The next version of Rubi will be extended to integrate expressions of the form

(g + h*x)^m * F^(f*(a + b*Log[c*(d + e*x)^n])^2)

Unfortunately Heyde's g(x) is not of that form.