Re: Random points in triangle
- To: mathgroup at smc.vnet.net
- Subject: [mg111617] Re: Random points in triangle
- From: Bill Rowe <readnews at sbcglobal.net>
- Date: Sat, 7 Aug 2010 06:21:53 -0400 (EDT)
On 8/6/10 at 6:55 AM, stevebg at 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.