Re: [R] generate random numbers subject to constraints

From: Charles C. Berry <cberry_at_tajo.ucsd.edu>
Date: Thu, 27 Mar 2008 08:17:22 -0700

On Thu, 27 Mar 2008, Robert A LaBudde wrote:

> At 05:06 PM 3/26/2008, Ted Harding wrote:
>> On 26-Mar-08 21:26:59, Ala' Jaouni wrote:
>>> X1,X2,X3,X4 should have independent distributions. They should be
>>> between 0 and 1 and all add up to 1. Is this still possible with
>>> Robert's method?
>>>
>>> Thanks
>>
>> I don't think so. A whileago you wrote
>> "The numbers should be uniformly distributed" (but in the
>> context of an example where you had 5 variable; now you
>> are back to 4 variables). Let's take the 4-case first.
>>
>> The two linear constraints confine the point (X1,X2,X3,X4)
>> to a triangular region within the 4-dimensional unit cube.
>> Say it has vertices A, B, C.
>> You could then start by generating points uniformly distributed
>> over a specific triangle in 2 dimentions, say the one with
>> vertices at A0=(0,0), B0=(0,1), C0=(1,0). This is easy.
>>
>> Then you need to find a linear transformation which will
>> map this triangle (A0,B0,C0) onto the triangle (A,B,C).
>> Then the points you have sampled in (A0,B0,C0) will map
>> into points which are uniformly distributed over the
>> triangle (A,B,C).
>>
>> More generally, you will be seeking to generate points
>> uniformly distributed over a simplex.
>>
>> For example, the case (your earlier post) of 5 points
>> with 2 linear constraints requires a tetrahedron with
>> vertices (A,B,C,D) in 5 dimensions whose coordinates you
>> will have to find. Then take an "easy" tetrahedron with
>> vertices (A0,B0,C0,D0) and sample uniformly within this.
>> Then find a linear mapping from (A0,B0,C0,D0) to (A,B,C,D)
>> and apply this to the sampled points.
>>
>> This raises a general question: Does anyone know of
>> an R function to sample uniformly in the interior
>> of a general (k-r)-dimensional simplex embedded in
>> k dimensions, with (k+1) given vertices?
>> <snip>
>
> The method of "rejection":
>
> 1. Generate numbers randomly in the hypercube.
> 2. Test to see if the point falls within the prescribed area.
> 3. Accept the point if it does.
> 4. Repeat if it doesn't.
>
> Efficiency depends upon the ratio of volumes involved.

The ratio is zero.

The subspace of the solution has lower dimension than the space you are sampling from.

So you will repeat '4' forever. (Up to machine accuracy, of course.)

And as I pointed out in my response to Ala' Jaouni, the 'solution' may lie in the null space. When it does not it will be a point, a line segment, a piece of a plane, or a 3 dimensional simplex.

Chuck

>
>
> ================================================================
> Robert A. LaBudde, PhD, PAS, Dpl. ACAFS e-mail: ral_at_lcfltd.com
> Least Cost Formulations, Ltd. URL: http://lcfltd.com/
> 824 Timberlake Drive Tel: 757-467-0954
> Virginia Beach, VA 23464-3239 Fax: 757-467-2947
>
> "Vere scire est per causas scire"
>
> ______________________________________________
> R-help_at_r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

Charles C. Berry                            (858) 534-2098
                                             Dept of Family/Preventive Medicine
E mailto:cberry_at_tajo.ucsd.edu	            UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901

R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Thu 27 Mar 2008 - 19:24:12 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Thu 27 Mar 2008 - 20:30:24 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive