Message Boards Message Boards

0
|
480 Views
|
16 Replies
|
3 Total Likes
View groups...
Share
Share this post:

Why getting complex number solution to expectation?

Posted 1 month ago

Can someone explain why I am getting a complex number for this computation of an expected value?

POSTED BY: Paul Studtmann
16 Replies

Again, this might not be helpful.

A similar bug is found here:

Integrate[Boole[R > 1/(1/S + 1/T) > P],
 {T, -1, 1}, {R, -1, 1}, {P, -1, 1}, {S, -1, 1}]

It returns an expression numerically equal to 2.88629 - 6.28319 I. Since Boole[] has the value either 0 or 1, the integral cannot be complex. And substituting NIntegrate[] for Integrate[] yields the real part 2.88629. The complex value returned by Integrate[] is the same as the value used in constructing the PDF of the transformed distribution for NExpectation[]/Expectation[] in the OP's problem. Interestingly it's not the same expression, and FullSimplify[] cannot transform their difference into zero. Nonetheless, they agree to 100 digits. (Probably, if I looked up the right PolyLog identity, I could do the simplification by hand. But I don't see the point.)

Tools for workarounds

If this is important, here are some ways to figure out the integral being used and how to fix it. The fixes are ad hoc and require some preliminary work (such as verifying the correct PDF).

The following return the integrals used for NExpectation[] and Expectation[] respectively, wrapped in Hold[] so they do no evaluate. If they look good, you could use ReleaseHold[] to evaluate them. Otherwise you would want to fix them before evaluating them.

integralNExp = Trace[
  NExpectation[R \[Conditioned] R > 1/(1/S + 1/T) > P,
   {T \[Distributed] UniformDistribution[{-1, 1}], 
    R \[Distributed] UniformDistribution[{-1, 1}], 
    P \[Distributed] UniformDistribution[{-1, 1}], 
    S \[Distributed] UniformDistribution[{-1, 1}]}],
  i_NIntegrate :> Return[Hold[i], Trace],
  TraceInternal -> True]

integralExp = Trace[
  Expectation[R \[Conditioned] R > 1/(1/S + 1/T) > P,
   {T \[Distributed] UniformDistribution[{-1, 1}], 
    R \[Distributed] UniformDistribution[{-1, 1}], 
    P \[Distributed] UniformDistribution[{-1, 1}], 
    S \[Distributed] UniformDistribution[{-1, 1}]}],
  i_Integrate :> Return[Hold[i], Trace],
  TraceInternal -> True]

In the OP's case there is a term -8 I Pi that should be removed, which I do by replacing the complex number -8 I by 0. That this is valid was verified by the Integrate[]/NIntegrate[] computations discussed at the beginning of this reply.

Trace[
 NExpectation[R \[Conditioned] R > 1/(1/S + 1/T) > P,
  {T \[Distributed] UniformDistribution[{-1, 1}], 
   R \[Distributed] UniformDistribution[{-1, 1}], 
   P \[Distributed] UniformDistribution[{-1, 1}], 
   S \[Distributed] UniformDistribution[{-1, 1}]}],
 i_NIntegrate :> 
  Return[Hold[i] /. -8 I -> 0 // Echo // ReleaseHold, Trace],
 TraceInternal -> True]
(*  0.5000000106328868`  *)

This yields a convenient (?), 2.6MB (!!!) expression. Even though Echo[] reveals we zeroed out all the complex terms, Integrate[] gives us another complex result, just like with the first example above.

Trace[
 Expectation[R \[Conditioned] R > 1/(1/S + 1/T) > P,
  {T \[Distributed] UniformDistribution[{-1, 1}], 
   R \[Distributed] UniformDistribution[{-1, 1}], 
   P \[Distributed] UniformDistribution[{-1, 1}], 
   S \[Distributed] UniformDistribution[{-1, 1}]}],
 i_Integrate :> 
  Return[Hold[i] /. -8 I -> 0 // Echo // ReleaseHold, Trace],
 TraceInternal -> True]

N[%, 16]
(*  0.500000000000000 + 2.176903850077749 I  *)

The real part agrees with the NExpectation[]/NIntegrate[] result, so maybe it's correct.

POSTED BY: Michael Rogers
Posted 26 days ago

My speculation about the appearance of complex numbers is that Mathematica uses the pdf of 1/(1/s+1/t) which contains logs:

dist = TransformedDistribution[1/(1/s + 1/t), {t \[Distributed] UniformDistribution[{-1, 1}], 
    s \[Distributed] UniformDistribution[{-1, 1}]}];
pdf = PDF[dist, st] // FullSimplify

pdf of 1/(1/s+1/t)

Because that random variable can take negative values and there are logs in the pdf, the internal code finds complex numbers even though those complex numbers will eventually "cancel out" outside of that internal code.

So a brute force approach to get the correct result is to do the numerical integeration:

numerator = NIntegrate[r (1/2) (1/2) pdf, {r, -1, 1}, {st, -1, r}, {p, -1, st}]
(* 0.0901967 *)
denominator = NIntegrate[(1/2) (1/2) pdf, {r, -1, 1}, {st, -1, r}, {p, -1, st}]
(* 0.180393 *)
conditionalMean = numerator/denominator
(* 0.5 *)

If one has a few minutes of patience, then Integrate also works:

numerator = Integrate[r (1/2) (1/2) pdf, {r, -1, 1}, {st, -1, r}, {p, -1, st}] // FullSimplify
(* 1/64 (3 + Log[16]) *)
denominator = Integrate[(1/2) (1/2) pdf, {r, -1, 1}, {st, -1, r}, {p, -1, st}] // FullSimplify
(* 1/32 (3 + Log[16]) *)
conditionalMean = numerator/denominator
(* 1/2 *)

Note: If Log[(-1 + st)/st] isn't split into two terms Log[-1 + st] - Log[st], then the imaginary numbers are not an issue. This condition happens when st < -1/2.

POSTED BY: Jim Baldwin
Posted 26 days ago

I am using 14.0 and am still getting weird results. I tried with Expectation rather than NExpectation. Using Expectation eliminates the 'invalid input' but still gives a complex number. When I compute:

Expectation[
 R \[Conditioned] (P < 1/(1 + 1/T) && 
    1/(1 + 1/T) < R), {T \[Distributed] UniformDistribution[{-1, 1}], 
  R \[Distributed] UniformDistribution[{-1, 1}], 
  P \[Distributed] UniformDistribution[{-1, 1}], 
  S \[Distributed] UniformDistribution[{-1, 1}]}]

I get a real number. But when I compute:

Expectation[
 R \[Conditioned] 
  R > 1/(1/S + 1/T) > P, {T \[Distributed] 
   UniformDistribution[{-1, 1}], 
  R \[Distributed] UniformDistribution[{-1, 1}], 
  P \[Distributed] UniformDistribution[{-1, 1}], 
  S \[Distributed] UniformDistribution[{-1, 1}]}]

I get a complex number. The only difference is that the first condition has 1/(1+1/T) where the second condition has 1/(1/S+1/T). So something seems amiss.

POSTED BY: Paul Studtmann
Posted 26 days ago

All your examples works very well with the right order on Mathematica 13.2:

enter image description here

POSTED BY: Denis Ivanov

I don't know if the following is any help:

This syntax seems valid (removing /T from first argument):

NExpectation[
 R \[Conditioned] R > 1/(1 + 1/S) > P,
 {T \[Distributed] UniformDistribution[{-1, 1}],
   R \[Distributed] UniformDistribution[{-1, 1}], 
  P \[Distributed] UniformDistribution[{-1, 1}], 
  S \[Distributed] UniformDistribution[{-1, 1}]}]

But the original throws an error. Maybe it's a bug?

Similarly this is accepted:

NExpectation[
 R \[Conditioned] R > T > S > P,
 {T \[Distributed] UniformDistribution[{-1, 1}], 
  R \[Distributed] UniformDistribution[{-1, 1}], 
  P \[Distributed] UniformDistribution[{-1, 1}], 
  S \[Distributed] UniformDistribution[{-1, 1}]}]

But this is not:

NExpectation[
 R \[Conditioned] R > S > P,
 {T \[Distributed] UniformDistribution[{-1, 1}], 
  R \[Distributed] UniformDistribution[{-1, 1}], 
  P \[Distributed] UniformDistribution[{-1, 1}], 
  S \[Distributed] UniformDistribution[{-1, 1}]}]

Perhaps there is a limitation on the number of variables that can be used in the first argument. Maybe the limitation is due to a bug?

I suggest reporting it to Wolfram. They may be able to help. If it is a bug, it needs to be fixed.

POSTED BY: Michael Rogers
Posted 26 days ago

Screenshot: enter image description here

POSTED BY: Denis Ivanov
Posted 26 days ago

What about 13.2?
Screenshot:

enter image description here

IMHO version 13.2 is the most stable, tested and balanced of the most recent.

POSTED BY: Denis Ivanov
Posted 26 days ago

Not only parentheses, but more correct order of inequality

POSTED BY: Denis Ivanov

I tried your modification in versions 12.0, 14.0, 14.1 and they all give the same complex number. Any idea why (a<b<c) and c>b>a should give different answers?

POSTED BY: Gianluca Gorni

Modifications are

R \[Conditioned] (P < 1/(1/S + 1/T) < R)

instead of

R\[Conditioned]R>1/(1/S+1/T)>P

I double checked this on version 13.2 and it works

POSTED BY: Gianluca Gorni

@Denis's modification was to put parentheses around the condition. (I think.)

However, I, too, still get the "invalid input" error, which seems more significant than the complex number. Basic GIGO: You can't expect the output to make sense if the input does not.

POSTED BY: Michael Rogers
Posted 26 days ago

Modifications are

R \[Conditioned] (P < 1/(1/S + 1/T) < R)

instead of

R\[Conditioned]R>1/(1/S+1/T)>P

I double checked this on version 13.2 and it works

POSTED BY: Denis Ivanov
Posted 26 days ago

Same here. I am not sure what the modifications are supposed to be. When I cut and paste the code, I get the same complex number.

POSTED BY: Paul Studtmann

What are the modifications? I get the same complex number also in the modified version. I tried with version 12, and I get the complex number but without error messages.

POSTED BY: Gianluca Gorni
Posted 27 days ago

Thank you!

POSTED BY: Paul Studtmann
Posted 27 days ago

After some modifications:

NExpectation[
 R \[Conditioned] (P < 1/(1/S + 1/T) < R), {T \[Distributed] 
   UniformDistribution[{-1, 1}], 
  R \[Distributed] UniformDistribution[{-1, 1}], 
  P \[Distributed] UniformDistribution[{-1, 1}], 
  S \[Distributed] UniformDistribution[{-1, 1}]}]

Out: 0.5
POSTED BY: Denis Ivanov
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract