I would like to ask for kind help with a following problem using TransformedDistribution in Mathematica. First here is the beginning of the code, that works:

n = 5;
k = 4;
f1 = OrderDistribution[{UniformDistribution[{0, 1}], n}, k]
f2 = OrderDistribution[{f1, n}, k]
w[v_?NumberQ] := NIntegrate[(v - b2)*PDF[f2, b2], {b2, 0, v}]
b[v_?NumberQ] := v - w[v]

The above code works so far, and getting any evaluated number from function b[v] works as well, for example when I write b[0.5] it gives me a number 0.49981412080261206, which is ok. I will not explain the context, but what I want to do is get a transformed distribution:

f3 := TransformedDistribution[b[v], v \[Distributed] f1]

So I want to transform each point in f1 by the function b[v]. Unfortunately this doesn't work, since when I try to plot pdf of f3 I get an empty graph, also checking the domain of f3 gives me [-inf,+inf], while clearly the domain here should be [0,1], since the basic function f1 is on [0,1], and the transformation is not decreasing each point by more than it's value. In short, there must be something wrong with the way I am using TransformedDistribution in this context. I am guessing the problem arises, because of how the integral limit depends on the parameter itself. I would very much appreciate if I could get any help on how to fix this issue, and how to make it work, so that I can get the transformed distribution.

  • $\begingroup$ try w[v_] := Integrate[(v - b2)*PDF[f2, b2], {b2, 0, v}, Assumptions -> {0 < v <= 1}] b[v_] := Simplify[v - w[v]]? $\endgroup$
    – kglr
    Jun 3, 2016 at 1:50
  • $\begingroup$ No this doesn't work, I get an error: NIntegrate::nlim: b2 = v is not a valid limit of integration. >>, when trying to plot PDF of f3 $\endgroup$
    – olga
    Jun 3, 2016 at 2:22
  • $\begingroup$ Integrate (not NIntegrate). And, you should clear older definitions before redefining w and b (i.e. use ClearAll[w,b] before w[v_] := .... $\endgroup$
    – kglr
    Jun 3, 2016 at 2:48
  • 1
    $\begingroup$ Rather than diving into providing lots of Mma code, why don't you explain what your problem is? It certainly appears 'unusual', taking an order statistic of an order statistic, and then a transformation of an expectation ... it would be helpful to see the problem expressed statistically. $\endgroup$
    – wolfies
    Jun 3, 2016 at 4:40

1 Answer 1

n = 5;
k = 4;

f1 = OrderDistribution[{UniformDistribution[{0, 1}], n}, k];

f2 = OrderDistribution[{f1, n}, k];

The definition of w[v] can be done exactly with Integrate rather than NIntegrate

w[v_] = Assuming[Element[v, Reals],
   Integrate[(v - b2)*PDF[f2, b2], {b2, 0, v}] // Simplify];

b[v_] = v - w[v] // Simplify

enter image description here

f3 = TransformedDistribution[b[v], v \[Distributed] f1];

The mean and standard deviation for f3 are

{m, s} = (#[f3] // N) & /@ {Mean, StandardDeviation}

(*  {0.628867, 0.140471}  *)

The PDF for f3 is zero outside of the interval

PDF[f3, x][[1, 1, -1]]

(*  0 < x < 51351184/66927861  *)

% // N

(*  0. < x < 0.767262  *)

Plot[Evaluate@PDF[f3, x], {x, 0, 51351184/66927861}, Frame -> True, 
 Mesh -> {Range[m - 4 s, m, s], None}, 
 MeshStyle -> {Red, AbsolutePointSize[4]}, Exclusions -> True]

enter image description here

NIntegrate[PDF[f3, x], {x, 0, 51351184/66927861 - 10^-75}, 
 WorkingPrecision -> 50, MaxRecursion -> 100]

(*  0.99999999999999999999999999999879953560281933889629  *)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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