Re: Random points in triangle
- To: mathgroup at smc.vnet.net
- Subject: [mg111635] Re: Random points in triangle
- From: "S. B. Gray" <stevebg at ROADRUNNER.COM>
- Date: Sun, 8 Aug 2010 07:21:59 -0400 (EDT)
- References: <i3jc3s$cuh$1@smc.vnet.net>
- Reply-to: stevebg at ROADRUNNER.COM
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 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. > >