Mathematica 9 is now available
Services & Resources / Wolfram Forums
-----
 /
MathGroup Archive
2004
*January
*February
*March
*April
*May
*June
*July
*August
*September
*October
*November
*December
*Archive Index
*Ask about this page
*Print this page
*Give us feedback
*Sign up for the Wolfram Insider

MathGroup Archive 2004

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

Search the Archive

Re: Re: Re: normal distribution random number generation

  • To: mathgroup at smc.vnet.net
  • Subject: [mg51323] Re: [mg51283] Re: [mg51265] Re: normal distribution random number generation
  • From: sean kim <sean_incali at yahoo.com>
  • Date: Thu, 14 Oct 2004 06:36:38 -0400 (EDT)
  • Reply-to: sean_incali at yahoo.com
  • Sender: owner-wri-mathgroup at wolfram.com

i just took a look at Part[plotvec, 24]

the respective values, for me, are, 

375 and 686.  

so i did a test with the following. 

in the fresh kernel, I get exactly the same thing as
you have. 

but I decided to evaluate the code

In[5]:=
SeedRandom[1111111 ]; 
vec1=Table[Random[Real,{0,1000}],{10^6}]; 
vec2=Map[If[#<44,1,0]&,vec1]; 
ones=Flatten[Position[vec2,1]]; 
plotvec=Map[Length,
Split[Sort[Drop[ones-RotateRight[ones,1],1]]]]; 
Part[plotvec, 24]
ListPlot[Log[plotvec],PlotRange->All] 
ListPlot[plotvec,PlotRange->All] 

Out[10]=
375

and then 

reevaluate the codes 

In[13]:= SeedRandom[5] 
In[14]:= Timing[First[Table[Random[],{10000}]]]
Out[14]= {0.01 Second,0.786599}
In[15]:= <<RandomReplacement` 
Timing[First[Table[Random[],{10000}]]] 
Out[16]= {0.01 Second,0.611425}

and I get the different value than before(one with the
fresh kernel)

i wonder what's going on here?

sean 


--- Andrzej Kozlowski <akoz at mimuw.edu.pl> wrote:

> I think I have now fixed the problem reported by
> Mark Fisher. The new  
> version of the RandomReplacement now works like
> this:
> 
> 
> SeedRandom[5]
> 
> 
> Timing[First[Table[Random[],{10000}]]]
> 
> 
> {0.01 Second,0.786599}
> 
> 
> <<RandomReplacement`
> 
> Timing[First[Table[Random[],{10000}]]]
> 
> 
> {0.02 Second,0.409222}
> 
> 
> As you can see there is a loss of speed involved.
> There may also be  
> still some as yet undiscovered problems. However,
> the new package is  
> available form my sites, though I still should write
> some documentation  
> of the latest changes (a few bugs in the previous
> version have also  
> been fixed).
> 
> Andrzej Kozlowski
> 
> Andrzej Kozlowski
> Chiba, Japan
> http://www.akikoz.net/~andrzej/
> http://www.mimuw.edu.pl/~akoz/
> 
> 
> 
> On 12 Oct 2004, at 14:57, Andrzej Kozlowski wrote:
> 
> > *This message was transferred with a trial version
> of CommuniGate(tm)  
> > Pro*
> > It's clearly because of packed arrays:
> >
> >
> > <<Developer`
> >
> >
> > PackedArrayQ[Table[Random[],{249}]]
> >
> >
> > False
> >
> >
> > PackedArrayQ[Table[Random[],{250}]]
> >
> >
> > True
> >
> >
> >
> > Obviously Mathematica uses a different method of
> generating random
> > packed arrays. There are various ways of getting
> round this (e.g.
> > constructing lists as joins of lists of length
> less than 250 etc) but
> > one will loose the benefits of packed arrays and
> with that, presumably,
> > a lot of speed. This seems to be something that
> only WRI can change.
> > Note however that not all situations in which one
> uses Random[] require
> > generating this type of lists.
> > Moreover, for problems involving the normal
> distribution it seems to me
> > that using the Marsaglia generator gives good
> results even with the
> > built in Random[]. But certainly we should hope
> that something will
> > finally be done do deal with this issue in version
> 6.
> >
> > Andrzej
> >
> >
> >
> >
> > On 11 Oct 2004, at 17:24, DrBob wrote:
> >
> >> He's right! How the devil does that happen?
> >>
> >> Here's a test, with Andrzej's package loaded in
> my Init file.
> >>
> >> First, with n=250:
> >>
> >> Quit
> >> SeedRandom[5]
> >> Table[Random[],{250}];
> >> Last@%
> >>
> >> 0.107874
> >>
> >> Unprotect[Random];
> >> Clear@Random
> >> SeedRandom[5]
> >> Table[Random[],{250}];
> >> Last@%
> >>
> >> 0.107874
> >>
> >> Now, with n=249:
> >>
> >> Quit
> >> SeedRandom[5]
> >> Table[Random[],{249}];
> >> Last@%
> >>
> >> 0.656099
> >>
> >> Unprotect[Random];
> >> Clear@Random
> >> SeedRandom[5]
> >> Table[Random[],{249}];
> >> Last@%
> >>
> >> 0.0708373
> >>
> >> Bobby
> >>
> >> On Mon, 11 Oct 2004 01:25:24 -0400 (EDT), Mark
> Fisher
> >> <mark at markfisher.net> wrote:
> >>
> >>> FYI: I've just a little testing and I find that
> Mathematica ignors
> >>> the user
> >>> defined rules for Random in Table[Random[],{n}]
> when n >= 250.
> >>>
> >>> Andrzej Kozlowski wrote:
> >>>> The problem is not actually with the way
> Mathematica's
> >>>> NormalDistribution but with the uniformly
> distributed Random[]
> >>>> function.
> >>>> NormalDistribution itself is generated (in the
> package
> >>>> Statistics`NormalDistribution`) by means of the
> very classical so
> >>>> called Boox-Muller method (although actually
> the Marsaglia variant
> >>>> below works better). You can try downolading my
> little
> >>>> RandomReplacement package form one of my web
> sites:
> >>>>
> >>>> http://www.akikoz.net/~andrzej//Mathematica/ 
> in Japan
> >>>>
> >>>> or
> >>>>
> >>>> http://www.mimuw.edu.pl/~akoz/Mathematica/   in
> Poland
> >>>>
> >>>> The package is based on the post of daniel
> Lichtblau, which explains
> >>>> the nature of the problem:
> >>>>
> >>>>
>
<http://library.wolfram.com/mathgroup/archive/2000/May/
> 
> >>>> msg00088.html>
> >>>>
> >>>>
> >>>>
> >>>> The package is named init.m and can be loaded
> automatically at the
> >>>> start of each Mathematica session.
> >>>> However, in fact if you are only concerned with
> normal distributions
> >>>> it
> >>>> may be enough to the following. First load the
> (unmodified)
> >>>> NormalDistribution package using
> >>>>
> >>>> <<Statistics`NormalDistribution`
> >>>>
> >>>> and then evaluate the code below. It will
> replace the Box-Muller
> >>>> generator by the Marsaglia one. I have found
> that this is usually
> >>>> enough. But if your problem still remains try
> the RandomReplacement
> >>>> package. (Or you can use both of course). Here
> is the Marsaglia code
> >>>> for normal distribution:
> >>>>
> >>>>
> >>>> Statistics`NormalDistribution`Private`normal=
> >>>>     
>
Compile[{{mu,_Real},{sigma,_Real},{q1,_Real},{q2,_Real}},
> >>>>        Module[{va=1.,vb,rad=2.,den=1.},
> >>>>         
> While[rad>=1.,va=2.*Random[]-1.;vb=2.*Random[]-1.;
> >>>>            rad=va*va+vb*vb];
> >>>>          den=Sqrt[-2.*(Log[rad]/rad)];
> >>>>          mu+sigma*va*den]];
> >>>>
> >>>>
> >>>>
> >>>> On 5 Oct 2004, at 17:36, Chris Kunkel wrote:
> >>>>
> >>>>
> >>>>> *This message was transferred with a trial
> version of
> >>>>> CommuniGate(tm)
> >>>>> Pro*
> >>>>> Hello,
> >>>>>
> >>>>> Recently I have encountered a problem in
> Mathematica's normal
> >>>>> distribution random number generator.  The
> problem arises when I
> >>>>> look
> >>>>> at the distribution of averages of a list of
> theses numbers.  That
> >>>>> is,
> >>>>> I generate 1000 numbers and take their
> average.  I do this a number
> >>>>> of
> >>>>> times and the plot a frequency distribution. 
> Consistently it seems
> >>>>> to
> >>>>> be skewed positive.  Specifically, the number
> of occurrences less
> >>>>> than
> >>>>> 3 standard deviations is consistent with the
> expected number, but
> >>>>> the
> >>>>> number greater than 3 is always larger (in
> some cases about twice
> >>>>> the
> >>>>> expected number).
> >>>>>
> >>>>> I was wondering if anyone else has noticed
> this and knows of a fix.
> >>>>> Also, does anyone have code for a good quality
> normal distribution
> >>>>> random number generator?
> >>>>>
> >>>>> Chris
> >>>>>
> >>>>
> >>>>
> >>>
> >>>
> >>>
> >>>
> >>
> >>
> >>
> >> -- 
> >> DrBob at bigfoot.com
> >> www.eclecticdreams.net
> >>
> >
> >
> 
> 





		
__________________________________
Do you Yahoo!?
Read only the mail you want - Yahoo! Mail SpamGuard.
http://promotions.yahoo.com/new_mail 


  • Prev by Date: Re: Re: Re: normal distribution random number generation
  • Next by Date: Re: Re: Re: normal distribution random number generation
  • Previous by thread: Re: Re: Re: normal distribution random number generation
  • Next by thread: Re: Re: Re: normal distribution random number generation