# TransformedDistribution using integration with parameter doesn't work

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.

• try w[v_] := Integrate[(v - b2)*PDF[f2, b2], {b2, 0, v}, Assumptions -> {0 < v <= 1}] b[v_] := Simplify[v - w[v]]? – kglr Jun 3 '16 at 1:50
• 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 – olga Jun 3 '16 at 2:22
• Integrate (not NIntegrate). And, you should clear older definitions before redefining w and b (i.e. use ClearAll[w,b] before w[v_] := .... – kglr Jun 3 '16 at 2:48
• 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. – wolfies Jun 3 '16 at 4:40

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 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}, Exclusions -> True] NIntegrate[PDF[f3, x], {x, 0, 51351184/66927861 - 10^-75},
WorkingPrecision -> 50, MaxRecursion -> 100]

(*  0.99999999999999999999999999999879953560281933889629  *)