[Date Index] [Thread Index] [Author Index]

Re: RandomArray from user defined distribution?

• To: mathgroup at smc.vnet.net
• Subject: [mg73390] Re: RandomArray from user defined distribution?
• From: Roland Franzius <roland.franzius at uos.de>
• Date: Thu, 15 Feb 2007 04:59:32 -0500 (EST)
• Organization: Universitaet Hannover
• References: <eqpdv5$eip$1@smc.vnet.net> <equndr$idu$1@smc.vnet.net>

roland franzius wrote:


>> rob wrote:
>>> I'd like to use the RandomArray to produce some data from 
>>> what I think is a Poisson distribution in t
>>> P[t] = a^2 t Exp[-a*t]  where a is mean, sigma.
>>>
>>> I see one can use RandomArray to produce sample data from a 
>>> lot of continuous distributions but the Poisson isn't among 
>>> them (it's only available in the discrete form).
>>>
>>> I've made a bunch of crippled attempts to force my P[t] to 
>>> put out examples but have failed. Any suggestions? Thanks.
>>>
>> 
>> Your density is a GammaDistribution[2,3], se e.g.
>> 
>> http://mathworld.wolfram.com/GammaDistribution.html

 >Sir, I email you directly because I can't seem to get my posts 
accepted >every time on mathgroup.
 >
 >Thank you for the suggestion on the Gamma function which is definitely 
 >what I was looking for and had not found in my searches.
 >
 >It appears I need GammaDistribution[2, mean/2] to get what I want.
 >
 >On your other suggestions, I don't follow. Please, what do you mean by >?
 >choose:=ip[Random[Real,{0,1}]]

Hi,

(*Perhaps this short example will enable you to apply the method*)

(*an exotic distibution density*)

density[x_] = (#/Integrate[#, {x, 0, ?}] &)[Sin[x^2]^2*Exp[-x^2/12]]

Plot[density[x], {x, 0, 12}, PlotRange -> All, PlotLabel -> "density"];

(*its cumulative distribution function not invertible by algebraic 
methods *)

distribution[x_] = Integrate[density[y], {y, 0, x}]

(*the x -> F(x) plot*)

Plot[distribution[x], {x,0, 8},
PlotRange -> {0, 1}, PlotLabel -> "distibution"]

(* a 1024 point table of points of its inverse y -> F^(-1)(y),
chopping removes small imaginary parts, the small additional increment 
table separates points in the first argument*)

table =
Chop@N@Table[
{distribution[x], x}, {x, 0, 12, 12/1024}] +
	Table[{x*10^-7, 0}, {x, 0, 1, 1/1024}];

(* ListPlot showing the inverted graph*)

ListPlot[table, PlotJoined -> True,
     PlotLabel -> "inverse table", PlotRange -> {0, 10}];

(*construct an interpolating pure function Ip from that table*)

Ip = Interpolation[table]

(*shows the same smooth graph *)

Plot[Ip[x], {x, 0.0,
    1.0}, PlotLabel -> "inverse interpolation", PlotRange -> {0, 10}];

(*define a random variable from the uniform distribution on (0, 1) using \
Random[]*)
choose := Ip[Random[]]
sample = Table[choose, {10000}]
sorted = Sort[sample];

(* ListPlots sorted and unsorted show a smooth approximation of the 
distribution and the described randomness of raw sample data*)

sorted = Sort[Table[choose, {10000}]];
ListPlot[sorted, PlotLabel -> "sorted sample", PlotRange -> {0, 10}];
ListPlot[sample, PlotLabel -> "unsorted sample", PlotRange -> {0, 10}]


Hope it helps

-- 

Roland Franzius