Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I'm using Mathematica 9.0 to calculate the probability according to the empirical distribution function (EDF) of some sample data. Afterwords this is included in a maximization stage so I define this probability as a function. I have to do apply this to groups of data in which each group has different dimensionality. Therefore I'd like to define a generic version of "probability given by an EDF" that I can apply to any of these groups.

Here there is a toy example for some 2-dimensional samples:

data = Transpose[{{1, 3, 4, 9, 8, 7, 8}, {2, 1, 1, 6, 7, 8, 9}}];
MatrixForm[data]

Its EDF is simply:

edf := EmpiricalDistribution[data];

Then by hand one can easily define its probability function:

 edfProbFunction[t1_, t2_] := NProbability[x1 <= t1 \[And] x2 <= t2, {x1, x2} \[Distributed] edf];

and compute the probability by just defining:

edfProb[w1_, w2_] := Evaluate[edfProbFunction[w1, w2]];

In this way, given a new point (2,4) from this distribution, it has probability 0.142857 given by:

edfProb[2,4]

My question is how to define a function like edfProbFunction for any dimension, not a fixed dimension (2 in the example)?

I tried to do it in different ways but didn't succeed. I lack background in Mathematica so these attempts may be nonesense. I summarize them anyway in case this can be of any help:

First naive attempt -- use a vector of input variables, straightforward

edfProbFunction[t__] := NProbability[x <= t, x \[Distributed] edf];
edfProb[w__] := Evaluate[edfProbFunction[w]];

Using this definition of edfProb together with the edfProbFunction defined in the toy example works, but not with this edfProbFunction here. This made me think that I had to somehow made explicit each of the individual predicates (the inequalities) in edfProbFunction

Second attempt -- use MakeBoxes

But a simple example shows that the expressions produced by this are not seen as variables in edfProbFunction:

xP /: MakeBoxes[xP[x___], form_] := RowBox[Riffle[Map[MakeBoxes[#, form] &, {x}], ","]]
x[1] = "" <> {"x", IntegerString[1]};
x[2] = "" <> {"x", IntegerString[2]};
varx = xP[x[1],x[2]] (*this produces a list x1,x2 *)
edfProbFunction[t1_,t2_] := NProbability[x1 <= t1 \[And] x2 <= t2, {varx} \[Distributed] edf];

Third attempt -- I tried to define a function that recursively creates the predicate, but it doesn't work neither:

table = Table[x[i] <= t[i], {i, 2}];
g[n_] := If[Length[n] > 1, n[[1]] \[And] g[Drop[n, 1]], If[Length[n] == 1, n[[1]], 0]]
predicate = g[table]
edfProbFunction[t_] := NProbability[predicate, {t[1], t[2]} \[Distributed] edf];
edfProb[w__] := Evaluate[edfProbFunction[w]];
edfProb[{2, 4}]

Any suggestions will be welcome, specially complete answers to my question.

share|improve this question

1 Answer 1

up vote 3 down vote accepted

The biggest issue in trying to put this together is that NProbability mixed with EmpiricalDistribution doesn't seem to like array indexed variables. Building up symbols programatically seems to fix the issue.

edfP[data_?MatrixQ][t__?NumericQ] /; 
  Length[{t}] == Length[data[[1]]] := Block[{vars, x},
  vars = Table[Symbol["x" <> ToString[i]], {i, Length[{t}]}]; 
  NProbability @@ {And @@ Thread[vars <= {t}], 
    vars \[Distributed] EmpiricalDistribution[data]}
  ]

First lets verify it works for your example...

edfProb[3, 4]
(* 0.285714 *)

edfP[data][3, 4]
(* 0.285714 *)

And now to generalize...

edfP[RandomVariate[NormalDistribution[], {10^4, 4}]][-1, 1, 2, 1]

(* 0.1061 *)

Note that if you are going to run this repeatedly for the same data you should evaluate the empirical distribution outside and pass that in rather than the data itself.

EDIT:

Now all that said, this seems like overkill to me. If I understand your question correctly why can't you just use CDF?

SeedRandom[1];
dat = RandomVariate[NormalDistribution[], {10^4, 4}];

edfP[dat][-1, 1, 2, 1]
(* 0.1143 *)

CDF[EmpiricalDistribution@dat, {-1, 1, 2, 1}]
(* 0.1143 *)
share|improve this answer
    
Both your approaches work really fine, they answer my question. So why can't I just use CDF? Well, I couldn't because I didn't even know this simple way, but now I can. Thanks a lot! –  p-d Mar 6 at 17:11

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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