Re: A conditional random number generation problem (please help me!)
- To: mathgroup at smc.vnet.net
- Subject: [mg65677] Re: A conditional random number generation problem (please help me!)
- From: dh <dh at metrohm.ch>
- Date: Fri, 14 Apr 2006 04:32:03 -0400 (EDT)
- References: <e1ijuv$3oh$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Hi,
if I understand correctly, you want a random number generator that gives
numbers in [a,b] with a given density P.
This is a bit tricky. First note that to make life easier we can assume
that the interval is [0,1]. This can then be mapped to [a,b] by:
a+x(b-a)
Further, Random[] gives random numbers in [0,1] uniformly distributed
with density c (==constant). We are looking for a function F[x] that
mappes the uniform distribution into one with the given density P[x]:
y=F[x]: [0,1] -> [0,1]
This must now be cast into mathematical notation. Assume we have a samll
intervall dx that contains dn random numbers. This will be mapped into
the interval:
dy== F'[x] dx
containing the same dn random numbers. Therefore,
dn== P[y] dy == P[F[x]] F'[x] dx == c dx
or
F'[x] == c/ P[ F[x] ]
That is the differential equation for F.
We also want that F mappes 0-> 0 and 1->1
Let's make an example:
we want a triangle distribution in [0,1]:
P[x_]:= 4 c If[x<0.5, x, 1-x]
the 4 c ensures that the numbers of random points stays constant. As
this function has a kink at 0.5, we integrate from 0..0.5 and get the
rest by symmetry. There is a further difficulty here, because P[0]==0,
and 1/P[0] is indefinite. We circumvent this problem by replacing 0 by a
very small number, say 10^-6.
Therefore, in [0,0.5] P[x]==x and c/P[F[x]]] == 0.25/F[x] and we can
solve the differential equation:
fun0 = F /. NDSolve[{F'[x] == 0.25/F[x], F[0] == 10^-6}, F,
{x, 0,0.5}][[1]];
fun[x_] := If[x < .5, fun0[x], 1 - fun0[1 - x]]
Plot[fun[x], {x, 0, 1}]
fun[x] will do the looked for mapping. We can test our result by simulation:
<< Graphics`Graphics`
d = Table[fun[Random[Real, {0, 1}]], {10000}];
Histogram[d]
Daniel
rjmachado3 wrote:
> I need to know the formula for the random function that return random
> numbers in a range of a and b integers [a,b] but that obey on a custom
> probability (possibly different!) for each integer number on this [a,b]
> range (of course the sum of all integer number probabilities are = 1!).
> Finally, what i want is the general function formula that simulate the
> random behavior (based on a custom probability value for each integer
> number between the [a,b] range. confuse? i hope not! please help me!!!!
>
> what i know so far is that the function formula for generating a "pure"
> random number between [a,b] range is:
>
> rand()*(b-a)+a
>
> where rand() return a random number between 0 and 1.
>