1
$\begingroup$

I am studying this piecewise function:

Piecewise[{{f1[x, h], x < a}, {f2[x, a, h], a <= x < b},     
{f3[x, a, b, h], b <= x < c}, {f4[x, a, b, c, h],  x >= c}}]

where

f1[x_,h_] := pdf1[x]/pdf[x] (h + (1 - h) 2 cdf1[x])/(h + (1 - h) 2 cdf[x])
f2[x_,a_,h_] := pdf1[x]/pdf[x] (h + (1 - h) 2 (1 - cdf1[x] + cdf1[a]))/
(h + (1 - h) 2 (1 - cdf[x] + cdf[a]))
f3[x_,a_,b_,h_] := pdf1[x]/pdf[x] (h + (1 - h) 2 (cdf1[x] - cdf1[b] + cdf1[a]))/
(h + (1 - h) 2 (cdf[x] - cdf[b] + cdf[a]))
f4[x_,a_,b_,c_] := pdf1[x]/pdf[x] (h + (1 - h) 2 (1 - cdf1[x] + cdf1[c] - cdf1[b] + 
cdf1[a]))/(h + (1 - h) 2 (1 - cdf[x] + cdf[c] - cdf[b] + cdf[a]))

where, in turn

pdf1[x_] := PDF[NormalDistribution[1, 1], x];
pdf[x_] := PDF[NormalDistribution[0, 1], x];
cdf1[x_] := CDF[NormalDistribution[1, 1], x];
cdf[x_] := CDF[NormalDistribution[0, 1], x];

I want the first and the third piece of the piecewise function to be smaller than a certain value $l$, say $l=0.48$; and the second and fourth piece of the function to be greater than the same $l$. The following code represents the intervals of $a$, $b$ and $c$ such that the four inequalities are satisfied:

{amin = NArgMin[{f2[a], f2[a] > l}, a], amax = NArgMax[{f1[a], f1[a] < l}, a], 
a = RandomReal[{amin, amax}]}
{bmin = a, bmax = NArgMax[{f3[b], f3[b] < l}, b], b = RandomReal[{bmin, bmax}]}
{cmin = b, cmax = NArgMax[{f3[c], f3[c] < l}, c], c = RandomReal[{cmin, cmax}]}

Note however that, for example, you have to give a random value to $a$ since the lower bound of the interval of the $b$'s depends on $a$.

I would like to plot in a 3D Graph all the possible combinations $(a,b,c)$ such that the system of inequalities is satisfied. Moreover, I would like to implement the Manipulate command so that I can study how the set of all possible combinations $(a,b,c)$ changes as $h$ varies between 0 and 1. If a 3D graph was not feasible, actually any graphical representation of the set of solutions of $(a,b,c)$ would be fine anyway. Thank you if you can help.

$\endgroup$
1
$\begingroup$

The requested plot can be obtained as follows. First, compute the limits on a, using the expressions in the question.

{amin = NArgMin[{f2[a], f2[a] > l}, a], amax = NArgMax[{f1[a], f1[a] < l}, a]}
(* {-0.233969, 0.544866} *)

Next, compute the corresponding limits in b.

btab = Interpolation@Table[{a, NArgMax[{f3[b], f3[b] < l}, b]}, {a, amin, amax, .02}];

and finally in c.

ctab = Interpolation[Flatten[Table[{{a, b}, Quiet@NArgMax[{f3[c], f3[c] < l, b < c},
    c]}, {a, amin, amax, .02}, {b, a, btab[a], .02}], 1], InterpolationOrder -> 1];

Since a < b < c, the plot is given by

Plot3D[{b, ctab[a, b]}, {a, amin, amax}, {b, a, btab[a]}, AxesLabel -> {a, b, c}, 
    LabelStyle -> {Black, Bold, Medium}, BoxRatios -> Automatic, ImageSize -> Large, 
    Mesh -> None, PlotStyle -> Opacity[.5]]

enter image description here

where the allowed values of c lie between the upper and lower curves (both of which have the same projection on the {a, b} plane). The limits on a, b, and c were obtained here.

Addendum

As noted by kgir in a comment above, RegionPlot3D also can be used once the logical structure and interpolation functions have been determined. RegionPlot3D is much slower here than Plot3D but produces a figure that some may prefer. The choice is a matter of taste.

Quiet@RegionPlot3D[amin < a < amax && a < b < btab[a] && b < c < ctab[a, b], 
    {a, amin, amax}, {b, amin, btab[amin]}, {c, amin, ctab[amin, btab[amin]]}, 
    AxesLabel -> {a, b, c}, LabelStyle -> {Black, Bold, Medium}, 
    BoxRatios -> Automatic, ImageSize -> Large, PlotPoints -> 500, Mesh -> None]

enter image description here

Second Addendum

As noted by the OP in a comment below, the computation of ctab given here often fails with the error

Interpolation::femimq The element mesh has insufficient quality of 0.`.

as h is increased. Basically, the triangular grid of Interpolation becomes tangled near the edge at which the two curves in the first plot come together. A simply way of eliminating the error is to replace 0.02 by 0.06 in the computation of ctab. However, doing so decreases the accuracy of the plot significantly for 1/3 < h < 1. A much more accurate but slower computation involves replacing the three lines of code immediately preceding the first plot by

btab = Interpolation@Table[{a, NArgMax[{f3[b], f3[b] < l}, b]}, {a, amin, amax}, 
    (amax - amin)/20}];
fc[a0_, b0_] := Quiet@NArgMax[{f3[c], f3[c] < l, b < c} /. {a -> a0, b -> b0}, c]
Plot3D[{b, fc[a, b]}, {a, amin, amax}, {b, a, btab[a]}, AxesLabel -> {a, b, c}, 
    LabelStyle -> {Black, Bold, Medium}, BoxRatios -> Automatic, ImageSize -> Large, 
    Mesh -> None, PlotStyle -> Opacity[.5], PlotPoints -> 10, MaxRecursion -> 0]

The resulting plot for h = 1/2 is

enter image description here

$\endgroup$
  • $\begingroup$ Thank you very much! I have a couple of questions: which value of $h$ did you use to get that plot? The limits of $a$, $b$ and $c$ depends on the $h$ you use. Is there a way to implement something like a Manipulate command so that I can study how the graph changes as $h$ varies? $\endgroup$ – Api Oct 5 '18 at 22:40
  • 1
    $\begingroup$ I used h = 0, but other values of h undoubtedly would work. Manipulate would be impractical in this case, because the computation requires several minutes, especially for the second plot. Please note that I added b < c to NArgMax[{f3[c], f3[c] < l, b < c}, c] from the other question, which I also updated, to treat unusual cases where f3 is not continuous. Incidentally, I soon will replace the two figures with ones of higher quality. $\endgroup$ – bbgodfrey Oct 5 '18 at 22:54
  • $\begingroup$ Ok, thank you for your clarification! I will just produce several graphs for different values of h and then try to come up with some nice way of presenting them. I will check for the higher quality figures then. Thank you for the time you are spending helping me. $\endgroup$ – Api Oct 6 '18 at 9:43
  • 1
    $\begingroup$ @Api The problem is with the ctab Interpolation grid near the edge where the two surfaces meet. A quick fix is to replace 0.02 by 0.06 in the definitions of btab and ctab, but the plot less accurate for larger h.. I would like to find a better solution, however. $\endgroup$ – bbgodfrey Oct 20 '18 at 2:10
  • 1
    $\begingroup$ @Api I have provided the code, but it is slow. $\endgroup$ – bbgodfrey Dec 9 '18 at 22:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.