Prev: ListPlot with PlotStyle->Dashed causes problems
Next: Nonlinear fitting question
From: S. B. Gray on 6 Aug 2010 06:55
I was looking for a simple way to place random points inside a triangle
with uniform distribution. Here's a good way:

newtri := Module[{x},
ptri = RandomReal[{-5, +5}, {3, 2}];
tredg = Subsets[ptri, {2}];
]
newpts[nump_] := Module[{wts},
inpoints = {};
Do [ wts = RandomReal[GammaDistribution[1, 2], 3];
wts = wts/Total[wts];
newin = Total[ptri*wts];
inpoints = Append[inpoints, newin], {nump}];
]
shotri := Module[{x},
Graphics[{Blue, Line[tredg], Red, Point[inpoints]}, ImageSize -> 500]
]

The same idea works for points in a tetrahedron; they will be uniformly
distributed if you use args such as GammaDistribution[.6,.1].

Steve Gray

From: Bill Rowe on 7 Aug 2010 06:21
On 8/6/10 at 6:55 AM, stevebg(a)ROADRUNNER.COM (S. B. Gray) wrote:

>I was looking for a simple way to place random points inside a
>triangle with uniform distribution. Here's a good way:

>newtri := Module[{x},
>ptri = RandomReal[{-5, +5}, {3, 2}]; tredg = Subsets[ptri, {2}];
>]
>newpts[nump_] := Module[{wts},
>inpoints = {}; Do [ wts = RandomReal[GammaDistribution[1, 2], 3];
>wts = wts/Total[wts]; newin = Total[ptri*wts]; inpoints =
>Append[inpoints, newin], {nump}];
>]
>shotri := Module[{x},
>Graphics[{Blue, Line[tredg], Red, Point[inpoints]}, ImageSize -> 500]
>]

>The same idea works for points in a tetrahedron; they will be
>uniformly distributed if you use args such as
>GammaDistribution[.6,.1].

It seems to me you are making something that should be
reasonably simple complex and doing so in a manner where you
aren't truly getting your stated desired result. Here, I take
your "random points inside a triangle with uniform distribution"
to mean points uniformly distributed over a specified triangle.

GammaDistribution[1, 2] is the same as ExponentialDistribution[2]

The sum of n independent items drawn from
ExponentialDistribution[b] will be distributed as
GammaDistribution[n, b]

The ratio x/(x+y) will have a beta distribution when x and y are
independent values drawn from gamma distributions.

A uniform distribution can be regarded as a special case of the
beta distribution with a = b = 1.

So, your variable newwin which is the product of two independent
things each coming from a beta distribution will certainly not
be uniformly distributed. However, the actual distribution for
newwin *might* be close enough to a uniform distribution for
your purpose. I haven't taken the time to work out the actual
distribution of newwin.

But rather than going through all of the above and trying to
work out the exact distribution for newwin, there is a simpler
approach which is guaranteed to give points uniformly
distributed over a given shape.

To select random points over some planar shape and ensure the
points are uniformly distributed over that shape, simply
generate a pair of uniform deviates using RandomReal[{-1,1}, 2].
That gives you uniformly distributed points over a square. Now
all you need do to get uniform points over your desired shape is
to drop points that lie outside your desired shape.

For simplicity, assume your desire is to have points distributed
uniformly over a unit circle. Then the values returned by:

Cases[RandomReal[{-1, 1},{1000,2}],_?(Norm@#<=1&)]

will be a set of points uniformly distributed over a circle with
radius 1

This can be readily see by doing

ListPlot[Cases[RandomReal[{-1, 1}, {1000, 2}], _?(Norm@# <= 1 &)],
AspectRatio -> 1]

The same idea can be extended for solids. For example,

Cases[RandomReal[{-1, 1},{1000,3}],_?(Norm@#<=1&)]

will be a set of points uniformly distributed through out the
volume of a sphere with radius 1.

Here, I chose a sphere and circle as examples because it is very
simple to determine whether a randomly selected point lies
inside or not. But the idea of doing this selection is valid no
matter what shape is desired.

Also, by using Cases, I am not controlling how many points I get
inside the desired area/volume. If you need to have a specific
number of points inside your desired shape, you could do
something like:

In[12]:= Table[pt = RandomReal[1, 2];
While[Norm[pt] > 1, pt = RandomReal[1, 2]]; pt, {5}]

Out[12]= {{0.337431, 0.0873031}, {0.25394, 0.927761}, {0.566057,
0.391587}, {0.497942, 0.441549}, {0.0179655, 0.854459}}

which gives you a list of 5 points inside or on the circle with
radius 1.

And of course it should go with out saying, larger or smaller
shapes can easily be obtained by multiplying your selected
points by the appropriate scale factor.

There is one other aspect of this approach that should be
mentioned. There are more random numbers being generated than
are actually used. If the desired shape is a very small fraction
of unit square, there will be a large number of the generated
points that are simply thrown away, causing a loss of
efficiency. That is, doing

Cases[RandomReal[{-10, 10},{100000,3}],_?(Norm@#<=1&)]

will generate approximately the same number of points in a
circle with radius 1 as

Cases[RandomReal[{-1, 1},{1000,3}],_?(Norm@#<=1&)]

at the cost of doing about 100 times as much work.

As long as the area/volume of the desired shape is reasonably
close to a square/cube then this shouldn't be much of a problem.
In fact, this method may out perform your approach since it is
takes less computation to generate a set of uniform deviates
that it does to generated deviates from other distributions.


From: Ray Koopman on 7 Aug 2010 06:22
On Aug 6, 3:55 am, "S. B. Gray" <stev...(a)ROADRUNNER.COM> wrote:
> I was looking for a simple way to place random points inside
> a triangle with uniform distribution. Here's a good way:
>
> newtri := Module[{x},
> ptri = RandomReal[{-5, +5}, {3, 2}];
> tredg = Subsets[ptri, {2}];
> ]
> newpts[nump_] := Module[{wts},
> inpoints = {};
> Do [ wts = RandomReal[GammaDistribution[1, 2], 3];
> wts = wts/Total[wts];
> newin = Total[ptri*wts];
> inpoints = Append[inpoints, newin], {nump}];
> ]
> shotri := Module[{x},
> Graphics[{Blue, Line[tredg], Red, Point[inpoints]}, ImageSize -> 500]
> ]
>
> The same idea works for points in a tetrahedron; they will be uniformly
> distributed if you use args such as GammaDistribution[.6,.1].

The scale parameter in the call to GammaDistribution doesn't matter,
because the scale gets divided out when you take wts/Total[wts], but
the shape parameter should always be 1. Values other than 1 will give
either too many or too few points near the vertices, according as the
shape parameter is less than or greater than 1.

From: S. B. Gray on 8 Aug 2010 07:21
Bill,
Thank you for your careful reply. I will consider what you say. In the
meantime, if you run the code I sent, you will see that the triangle is
filled with a pretty good uniform distribution. Visual uniformity is
really all that I need, not any statistical tests.
Another factor is that I have convex 3D polyhedra to be filled
"uniformly." It turns out that the Gamma distribution does a decent job
with a good setting of the first parameter; the more points there are in
the convex polyhedron, the lower that parameter should be for good
uniformity. (There is one weight, 0<wts<1, per vertex and the resulting
point is always inside the polyhedron.)
The purpose for the random internal points (and random convex
polyhedra) is to test a conjecture that I thought of. The conjecture has
to do with the angles between the n rays from a random point inside the
convex polyhedron to its n vertices. There are n(n-1)/2 angles between
pairs of rays and the conjecture is that at least n-1 of them must be
obtuse. This may be easy or difficult to prove - I don't know, and if
you have any ideas I'd be glad to hear them.

Steve Gray

On 8/7/2010 3:21 AM, Bill Rowe wrote:
> On 8/6/10 at 6:55 AM, stevebg(a)ROADRUNNER.COM (S. B. Gray) wrote:
>
>> I was looking for a simple way to place random points inside a
>> triangle with uniform distribution. Here's a good way:
>
>> newtri := Module[{x},
>> ptri = RandomReal[{-5, +5}, {3, 2}]; tredg = Subsets[ptri, {2}];
>> ]
>> newpts[nump_] := Module[{wts},
>> inpoints = {}; Do [ wts = RandomReal[GammaDistribution[1, 2], 3];
>> wts = wts/Total[wts]; newin = Total[ptri*wts]; inpoints =
>> Append[inpoints, newin], {nump}];
>> ]
>> shotri := Module[{x},
>> Graphics[{Blue, Line[tredg], Red, Point[inpoints]}, ImageSize -> 500]
>> ]
>
>> The same idea works for points in a tetrahedron; they will be
>> uniformly distributed if you use args such as
>> GammaDistribution[.6,.1].
>
> It seems to me you are making something that should be
> reasonably simple complex and doing so in a manner where you
> aren't truly getting your stated desired result. Here, I take
> your "random points inside a triangle with uniform distribution"
> to mean points uniformly distributed over a specified triangle.
>
> GammaDistribution[1, 2] is the same as ExponentialDistribution[2]
>
> The sum of n independent items drawn from
> ExponentialDistribution[b] will be distributed as
> GammaDistribution[n, b]
>
> The ratio x/(x+y) will have a beta distribution when x and y are
> independent values drawn from gamma distributions.
>
> A uniform distribution can be regarded as a special case of the
> beta distribution with a = b = 1.
>
> So, your variable newwin which is the product of two independent
> things each coming from a beta distribution will certainly not
> be uniformly distributed. However, the actual distribution for
> newwin *might* be close enough to a uniform distribution for
> your purpose. I haven't taken the time to work out the actual
> distribution of newwin.
>
> But rather than going through all of the above and trying to
> work out the exact distribution for newwin, there is a simpler
> approach which is guaranteed to give points uniformly
> distributed over a given shape.
>
> To select random points over some planar shape and ensure the
> points are uniformly distributed over that shape, simply
> generate a pair of uniform deviates using RandomReal[{-1,1}, 2].
> That gives you uniformly distributed points over a square. Now
> all you need do to get uniform points over your desired shape is
> to drop points that lie outside your desired shape.
>
> For simplicity, assume your desire is to have points distributed
> uniformly over a unit circle. Then the values returned by:
>
> Cases[RandomReal[{-1, 1},{1000,2}],_?(Norm@#<=1&)]
>
> will be a set of points uniformly distributed over a circle with
> radius 1
>
> This can be readily see by doing
>
> ListPlot[Cases[RandomReal[{-1, 1}, {1000, 2}], _?(Norm@#<= 1&)],
> AspectRatio -> 1]
>
> The same idea can be extended for solids. For example,
>
> Cases[RandomReal[{-1, 1},{1000,3}],_?(Norm@#<=1&)]
>
> will be a set of points uniformly distributed through out the
> volume of a sphere with radius 1.
>
> Here, I chose a sphere and circle as examples because it is very
> simple to determine whether a randomly selected point lies
> inside or not. But the idea of doing this selection is valid no
> matter what shape is desired.
>
> Also, by using Cases, I am not controlling how many points I get
> inside the desired area/volume. If you need to have a specific
> number of points inside your desired shape, you could do
> something like:
>
> In[12]:= Table[pt = RandomReal[1, 2];
> While[Norm[pt]> 1, pt = RandomReal[1, 2]]; pt, {5}]
>
> Out[12]= {{0.337431, 0.0873031}, {0.25394, 0.927761}, {0.566057,
> 0.391587}, {0.497942, 0.441549}, {0.0179655, 0.854459}}
>
> which gives you a list of 5 points inside or on the circle with
> radius 1.
>
> And of course it should go with out saying, larger or smaller
> shapes can easily be obtained by multiplying your selected
> points by the appropriate scale factor.
>
> There is one other aspect of this approach that should be
> mentioned. There are more random numbers being generated than
> are actually used. If the desired shape is a very small fraction
> of unit square, there will be a large number of the generated
> points that are simply thrown away, causing a loss of
> efficiency. That is, doing
>
> Cases[RandomReal[{-10, 10},{100000,3}],_?(Norm@#<=1&)]
>
> will generate approximately the same number of points in a
> circle with radius 1 as
>
> Cases[RandomReal[{-1, 1},{1000,3}],_?(Norm@#<=1&)]
>
> at the cost of doing about 100 times as much work.
>
> As long as the area/volume of the desired shape is reasonably
> close to a square/cube then this shouldn't be much of a problem.
> In fact, this method may out perform your approach since it is
> takes less computation to generate a set of uniform deviates
> that it does to generated deviates from other distributions.
>
>


From: Ray Koopman on 8 Aug 2010 07:23
On Aug 7, 3:21 am, Bill Rowe <readn...(a)sbcglobal.net> wrote:
> On 8/6/10 at 6:55 AM, stev...(a)ROADRUNNER.COM (S. B. Gray) wrote:
>
>> I was looking for a simple way to place random points inside a
>> triangle with uniform distribution. Here's a good way:
>> newtri := Module[{x},
>> ptri = RandomReal[{-5, +5}, {3, 2}]; tredg = Subsets[ptri, {2}];
>> ]
>> newpts[nump_] := Module[{wts},
>> inpoints = {}; Do [ wts = RandomReal[GammaDistribution[1, 2], 3];
>> wts = wts/Total[wts]; newin = Total[ptri*wts]; inpoints =
>> Append[inpoints, newin], {nump}];
>> ]
>> shotri := Module[{x},
>> Graphics[{Blue, Line[tredg], Red, Point[inpoints]}, ImageSize -> 500]
>> ]
>> The same idea works for points in a tetrahedron; they will be
>> uniformly distributed if you use args such as
>> GammaDistribution[.6,.1].
>
> It seems to me you are making something that should be
> reasonably simple complex and doing so in a manner where you
> aren't truly getting your stated desired result. Here, I take
> your "random points inside a triangle with uniform distribution"
> to mean points uniformly distributed over a specified triangle.
>
> GammaDistribution[1, 2] is the same as ExponentialDistribution[2]
>
> The sum of n independent items drawn from
> ExponentialDistribution[b] will be distributed as
> GammaDistribution[n, b]
>
> The ratio x/(x+y) will have a beta distribution when x and y are
> independent values drawn from gamma distributions.
>
> A uniform distribution can be regarded as a special case of the
> beta distribution with a = b = 1.
>
> So, your variable newwin which is the product of two independent
> things each coming from a beta distribution will certainly not
> be uniformly distributed.

Try this a few times. It looks pretty uniform to me.

p = RandomReal[NormalDistribution[0,1],{3,2}];
r = #.p / Total[#,{2}] & @
RandomReal[ExponentialDistribution[1],{5000,3}];
ListPlot[{p,r},PlotRange->All,AspectRatio->1,Frame->True,Axes->None,
PlotStyle->{{Red,PointSize[.02]},{Black,PointSize[.005]}}]

> However, the actual distribution for
> newwin *might* be close enough to a uniform distribution for
> your purpose. I haven't taken the time to work out the actual
> distribution of newwin.
>
> But rather than going through all of the above and trying to
> work out the exact distribution for newwin, there is a simpler
> approach which is guaranteed to give points uniformly
> distributed over a given shape.
>
> To select random points over some planar shape and ensure the
> points are uniformly distributed over that shape, simply
> generate a pair of uniform deviates using RandomReal[{-1,1}, 2].
> That gives you uniformly distributed points over a square. Now
> all you need do to get uniform points over your desired shape is
> to drop points that lie outside your desired shape.
>
> For simplicity, assume your desire is to have points distributed
> uniformly over a unit circle. Then the values returned by:
>
> Cases[RandomReal[{-1, 1},{1000,2}],_?(Norm@#<=1&)]
>
> will be a set of points uniformly distributed over a circle with
> radius 1
>
> This can be readily see by doing
>
> ListPlot[Cases[RandomReal[{-1, 1}, {1000, 2}], _?(Norm@# <= 1 &)],
> AspectRatio -> 1]
>
> The same idea can be extended for solids. For example,
>
> Cases[RandomReal[{-1, 1},{1000,3}],_?(Norm@#<=1&)]
>
> will be a set of points uniformly distributed through out the
> volume of a sphere with radius 1.
>
> Here, I chose a sphere and circle as examples because it is very
> simple to determine whether a randomly selected point lies
> inside or not. But the idea of doing this selection is valid no
> matter what shape is desired.
>
> Also, by using Cases, I am not controlling how many points I get
> inside the desired area/volume. If you need to have a specific
> number of points inside your desired shape, you could do
> something like:
>
> In[12]:= Table[pt = RandomReal[1, 2];
> While[Norm[pt] > 1, pt = RandomReal[1, 2]]; pt, {5}]
>
> Out[12]= {{0.337431, 0.0873031}, {0.25394, 0.927761}, {0.566057,
> 0.391587}, {0.497942, 0.441549}, {0.0179655, 0.854459}}
>
> which gives you a list of 5 points inside or on the circle with
> radius 1.
>
> And of course it should go with out saying, larger or smaller
> shapes can easily be obtained by multiplying your selected
> points by the appropriate scale factor.
>
> There is one other aspect of this approach that should be
> mentioned. There are more random numbers being generated than
> are actually used. If the desired shape is a very small fraction
> of unit square, there will be a large number of the generated
> points that are simply thrown away, causing a loss of
> efficiency. That is, doing
>
> Cases[RandomReal[{-10, 10},{100000,3}],_?(Norm@#<=1&)]
>
> will generate approximately the same number of points in a
> circle with radius 1 as
>
> Cases[RandomReal[{-1, 1},{1000,3}],_?(Norm@#<=1&)]
>
> at the cost of doing about 100 times as much work.
>
> As long as the area/volume of the desired shape is reasonably
> close to a square/cube then this shouldn't be much of a problem.
> In fact, this method may out perform your approach since it is
> takes less computation to generate a set of uniform deviates
> that it does to generated deviates from other distributions.


 |  Next  |  Last
Pages: 1 2
Prev: ListPlot with PlotStyle->Dashed causes problems
Next: Nonlinear fitting question