# Produce a probability vector with uniform distribution [duplicate]

A probability vector is a vector $\mathbf{p}=(p_1,p_2,\ldots,p_n)$ with the following properties:

• $\sum_{i=1}^{n}p_i=1$

• $p_i\geq 0$

Assume that the set of all possible probability vectors is $\mathcal{P}$. I want to uniformly pick a $\mathbf{p}$ from $\mathcal{P}$. How can I do that in Mathematica?

• This is a site for the computing software Mathematica, not for MATLAB. – march Jul 12 '17 at 19:46

I don't understand why you use the term "uniform distribution" or what you mean by "pick uniformly a $p$ from $P$". None of the marginal distributions are uniform (except in maybe a degenerate case). And in terms of notation you're using the same notation for the parameters of the distribution and realizations from that distribution.

I assume you do want to pick independent random samples from a distribution with such restrictions that you list. One such distribution is the Dirichlet distribution (a multivariate beta distribution):

p1 = 0.3;
p2 = 0.2;
p3 = 1 - p1 - p2;
nSamples = 10;
x = RandomVariate[DirichletDistribution[{p1, p2, p3}], nSamples];
(* {{0.00646779, 0.000392286}, {0.726117, 0.0721349}, {0.476097,
0.523398}, {3.77399*10^-6, 0.000228199}, {0.142504,
0.402989}, {0.839344, 0.146529}, {0.987459, 0.0123547}, {0.096979,
0.239194}, {0.628002, 0.0148397}, {0.978951, 0.00772091}} *)


The above set of random samples have only two values where it's understood that the third value is 1 minus the sum of the other two. To include the third value something like the following could be used:

x = Table[Flatten[{x[[i]], 1 - Total[x[[i]]]}], {i, Length[x]}]
(*{{0.00646779, 0.000392286, 0.99314}, {0.726117, 0.0721349,
0.201748}, {0.476097, 0.523398, 0.000504177}, {3.77399*10^-6,
0.000228199, 0.999768}, {0.142504, 0.402989, 0.454507}, {0.839344,
0.146529, 0.0141275}, {0.987459, 0.0123547, 0.000186575}, {0.096979,
0.239194, 0.663827}, {0.628002, 0.0148397, 0.357158}, {0.978951,
0.00772091, 0.013328}} *)


If by "uniformly" you mean that the joint density function is a constant positive number (where it's not zero), then you would use

nSamples = 10;
n = 3;
x = RandomVariate[DirichletDistribution[Table[1, {i, n}]], nSamples];
x = Table[Flatten[{x[[i]], 1 - Total[x[[i]]]}], {i, Length[x]}]

rp[dim_, samples_:1]:=Append[#, 1-Total @ #] & /@ RandomPoint[Simplex[dim], samples];

{#, rp @ #} & /@ Range[5]


{#, rp[#, 3]} &/@ Range[5]


rp[30]


{{0.0547105,0.0234894,0.039695,0.0159973,0.0596655,0.0657402,0.00261608,0.012769,0.0922238,0.0333687,0.0220436,0.00159148,0.00332486,0.00770202,0.0486021,0.0147596,0.00980909,0.0386692,0.00453807,0.0492418,0.181455,0.0171443,0.0155691,0.00373707,0.00245814,0.0183238,0.00451847,0.0229032,0.0251677,0.072695,0.0354716}}

The default Simplex (as used by @kglr in his answer) includes the origin. So, for your question, we could instead use:

n = 4;
pts = RandomPoint[Simplex[IdentityMatrix[n]], 5]


{{0.146941, 0.172728, 0.377783, 0.302548}, {0.0262286, 0.559292, 0.367651, 0.046828}, {0.697806, 0.0495862, 0.020688, 0.231919}, {0.266115, 0.295599, 0.0533206, 0.384965}, {0.708444, 0.188191, 0.0475761, 0.0557884}}

Check:

Total[pts, {2}]


{1., 1., 1., 1., 1.}

• For whatever it's worth, your method is a whole lot faster for large numbers of samples (say, 5,000 or more) than using the DirichletDistribution. – JimB Jul 12 '17 at 22:01
• @JimBaldwin With the introduction of RandomPoint, a total probability newbie like me can produce an answer. – Carl Woll Jul 12 '17 at 23:14