[Numpy-discussion] Re: A random.normal function with stdev as array

[Robert Kern]

Eric Emsellem wrote:
> Hi,
> 
> I am trying to optimize a code where I derive random numbers many times
> and having an array of values for the stdev parameter.
> 
> I wish to have an efficient way of doing something like:
> ##################
> stdev = array([1.1,1.2,1.0,2.2])
> result = numpy.zeros(stdev.shape, Float)
> for i in range(len(stdev)) :
>   result[i] = numpy.random.normal(0, stdev[i])
> ##################

You can use the fact that the standard deviation of a normal distribution is a
scale parameter. You can get random normal deviates of varying standard
deviation by multiplying a standard normal deviate by the desired standard
deviation (how's that for confusing terminology, eh?).

  result = numpy.random.standard_normal(stdev.shape) * stdev

> In my case,  stdev can in fact be an array of a few millions floats...
> so I really need to optimize things.
> 
> Any hint on how to code this efficiently ?
> 
> And in general, where could I find tips for optimizing a code where I
> unfortunately have too many loops such as "for i in range(Nbody) : "
> with Nbody being > 10^6 ?

Tim Hochberg recently made this list:

"""
0. Think about your algorithm.
1. Vectorize your inner loop.
2. Eliminate temporaries
3. Ask for help
4. Recode in C.
5. Accept that your code will never be fast.

Step zero should probably be repeated after every other step ;)
"""

That's probably the best general advice. To get better advice, we would need to
know the specifics of the problem.

-- 
Robert Kern
robert.kern at gmail.com

"I have come to believe that the whole world is an enigma, a harmless enigma
 that is made terrible by our own mad attempt to interpret it as though it had
 an underlying truth."
  -- Umberto Eco