Sample uniformly from a simplex

Andreas

unread,

Dec 12, 2008, 6:55:52 AM12/12/08

to

I need to develop Mathematica code to sample uniformly from a unit n-dimensional simplex.

I came across a description of the problem at: http://geomblog.blogspot.com/2005/10/sampling-from-simplex.html

Specifically, I would like a uniform sample from the set

X = { (x1, x2, ..., xD) | 0 <= xi <= 1, x1 + x2 + ... + xD = 1}.

D is the dimension of the simplex.

So, the coordinates of any point on the simplex would sum to 1 and I need to sample points on the simplex.

geomblog's solution suggested:

Generating IID random samples from an exponential distribution by sampling X from [0,1] uniformly, and returning -log(X)).

Take n samples, then normalize.

This should result in a list of numbers which is a uniform sample from the simplex.

I've searched extensively for a Mathematica implementation of something like this, to no avail.

I keep trying different things but haven't made much headway.

Any suggestions for how to develop this (or an equivelant) in Mathematica much appreciated

A

Andreas

unread,

Dec 14, 2008, 7:33:34 AM12/14/08

to

I've kept thinking about this and came up with the following.

sampleSimplex[sampleSize_, dimensions_] :=
Module[{sample, sampleOutput},
sample = RandomReal[1, {sampleSize, dimensions}];
sampleOutput =
Table[1 / Total[sample[[i, All]]] *sample[[i, All]], {i, 1,
sampleSize}];
sampleOutput
]

example:

sampleSimplex[4, 3]

{{0.163972, 0.445817, 0.390211}, {0.470828, 0.291766,
0.237406}, {0.56663, 0.3813, 0.0520696}, {0.182715, 0.600876, 0.216409}}

This doesn't draw from an exponential distribution, but would this work?

Also, if it does work can I do it a better (faster, more elegant) way?

Mark Fisher

unread,

Dec 14, 2008, 7:37:23 AM12/14/08

to

On Dec 12, 6:55 am, Andreas <aa...@ix.netcom.com> wrote:
> I need to develop Mathematica code to sample uniformly from a unit n-dime=
nsional simplex.
>
> I came across a description of the problem at:http://geomblog.blogspot.co=

m/2005/10/sampling-from-simplex.html
>
> Specifically, I would like a uniform sample from the set
>
> X = { (x1, x2, ..., xD) | 0 <= xi <= 1, x1 + x2 + ... + xD = 1}.
>
> D is the dimension of the simplex.
>

> So, the coordinates of any point on the simplex would sum to 1 and I need=

to sample points on the simplex.
>
> geomblog's solution suggested:
>

> Generating IID random samples from an exponential distribution by samplin=

g X from [0,1] uniformly, and returning -log(X)).
>
> Take n samples, then normalize.
>

> This should result in a list of numbers which is a uniform sample from th=
e simplex.
>
> I've searched extensively for a Mathematica implementation of something l=

ike this, to no avail.
>
> I keep trying different things but haven't made much headway.
>

> Any suggestions for how to develop this (or an equivelant) in Mathematica=
much appreciated
>
> A

Hi. Try this (using Version 6 or 7):

ranSimp[d_] := (#/Plus @@ #) &[-Log[RandomReal[1, d]]]
ranSimp[d_, n_] := (#/Plus @@ #) & /@ (-Log[RandomReal[1, {n, d}]])

d is the "dimension" of the simplex and n is the number of draws.

--Mark

Asim

unread,

Dec 14, 2008, 7:37:58 AM12/14/08

to

On Dec 12, 6:55 am, Andreas <aa...@ix.netcom.com> wrote:

> I need to develop Mathematica code to sample uniformly from a unit n-dime=
nsional simplex.
>
> I came across a description of the problem at:http://geomblog.blogspot.co=