[Numpy-discussion] numarray: RandomArray2

[Todd Miller]

Jochen Küpper wrote:

>On Mon, 19 Aug 2002 09:27:35 -0400 Todd Miller wrote:
>
>Todd> Jochen Küpper wrote:
>
>>>Also note the "strange" docstring of range:
>>>,----
>>>| range(...)
>>>|     range([start,] stop[, step]) -> list of integers
>>>`----
>>>pointing straight toward its behavior.
>>>
>
>Ok, reworked the docstrings of module RandomArray2 a bit in order to
>(hopefully) clarify these issues. I have attached a diff, if it is ok
>I'll check it in.
>
Please do.

Todd

>
>>>So well, maybe someone with insight into RNG's can comment on the,
>>>mirroring issue?
>>>
>
>Anybody?
>
>Todd> It appears to me that the parameter order of randint again
>Todd> emulates "range". The fact that random_integers is not
>Todd> consistent with randint seem OK to me because random_integers
>Todd> appears to have been written expressly to tailor the calling
>Todd> sequence of randint.
>
>Put a comment about that in the docstring of random_integers.
>
>Greetings,
>Jochen
>
>
>------------------------------------------------------------------------
>
>? Packages/RandomArray2/build
>? Packages/RandomArray2/Lib/ChangeLog
>Index: Packages/RandomArray2/Lib/RandomArray2.py
>===================================================================
>RCS file: /cvsroot/numpy/numarray/Packages/RandomArray2/Lib/RandomArray2.py,v
>retrieving revision 1.3
>diff -u -r1.3 RandomArray2.py
>--- Packages/RandomArray2/Lib/RandomArray2.py	16 Aug 2002 10:44:21 -0000	1.3
>+++ Packages/RandomArray2/Lib/RandomArray2.py	22 Aug 2002 22:11:16 -0000
>@@ -5,6 +5,11 @@
> import math
> from types import *
> 
>+__doc__ = """Random number array generators.
>+
>+This module provides functions to generate arrays of random numbers.
>+"""
>+
> # Extended RandomArray to provide more distributions:
> # normal, beta, chi square, F, multivariate normal,
> # exponential, binomial, multinomial
>@@ -13,8 +18,8 @@
> ArgumentError = "ArgumentError"
> 
> def seed(x=0,y=0):
>-    """seed(x, y), set the seed using the integers x, y; 
>-    Set a random one from clock if  y == 0
>+    """Set the RNG seed using the integers x, y; 
>+    If |y| == 0 set a random one from clock.
>     """
>     if type (x) != IntType or type (y) != IntType :
>         raise ArgumentError, "seed requires integer arguments."
>@@ -30,14 +35,14 @@
> seed()
> 
> def get_seed():
>-    "Return the current seed pair"
>+    """Return the current seed pair"""
>     return ranlib.get_seeds()
> 
> def _build_random_array(fun, args, shape=[]):
>-# Build an array by applying function fun to
>-# the arguments in args, creating an array with
>-# the specified shape.
>-# Allows an integer shape n as a shorthand for (n,).
>+    """Build an array by applying function |fun| to the arguments in
>+    |args|, creating an array with the specified shape.
>+    Allows an integer shape n as a shorthand for (n,).
>+    """
>     if isinstance(shape, IntType): 
>         shape = [shape]
>     if len(shape) != 0:
>@@ -51,20 +56,28 @@
>         return s[0]
> 
> def random(shape=[]):
>-    "random(n) or random([n, m, ...]) returns array of random numbers"
>+    """random(n) or random([n, m, ...])
>+
>+    Returns array of random numbers in the range 0.0 to 1.0.
>+    """
>     return _build_random_array(ranlib.sample, (), shape)
> 
> def uniform(minimum, maximum, shape=[]):
>-    """uniform(minimum, maximum, shape=[]) returns array of given shape of random reals 
>-    in given range"""
>+    """uniform([minimum,], maximum[, shape])
>+
>+    Return array with shape |shape| of random Floats with all values
>+    > minimum and < maximum.    
>+    """
>     return minimum + (maximum-minimum)*random(shape)
> 
> def randint(minimum, maximum=None, shape=[]):
>-    """randint(min, max, shape=[]) = random integers >=min, < max
>-    If max not given, random integers >= 0, <min"""
>+    """randint([minimum,] maximum[, shape])
>+
>+    Return random integers >= mininimum, < maximum,
>+    if maximum not given, random integers >= 0, < minimum.
>+    """
>     if maximum is None:
>-        maximum = minimum
>-        minimum = 0
>+        maximum, minimum = minimum, 0
>     a = Numeric.floor(uniform(minimum, maximum, shape))
>     if isinstance(a, Numeric.ArrayType):
>         return a.astype(Numeric.Int32)
>@@ -72,79 +85,94 @@
>         return int(a)
>      
> def random_integers(maximum, minimum=1, shape=[]):
>-    """random_integers(max, min=1, shape=[]) = random integers in range min-max inclusive"""
>+    """Return array of random integers in interval [minimum, maximum]
>+
>+    Note that this function has reversed arguments. It is simply a
>+    redirection through randint, and use of that function (randint) is
>+    suggested.
>+    """
>     return randint(minimum, maximum+1, shape) 
>      
> def permutation(n):
>-    "permutation(n) = a permutation of indices range(n)"
>+    """A permutation of indices range(n)"""
>     return Numeric.argsort(random(n))
> 
> def standard_normal(shape=[]):
>-    """standard_normal(n) or standard_normal([n, m, ...]) returns array of
>-           random numbers normally distributed with mean 0 and standard
>-           deviation 1"""
>+    """standard_normal(n) or standard_normal([n, m, ...])
>+
>+    Returns array of random numbers normally distributed with mean 0
>+    and standard deviation 1.
>+    """
>     return _build_random_array(ranlib.standard_normal, (), shape)
> 
> def normal(mean, std, shape=[]):
>-        """normal(mean, std, n) or normal(mean, std, [n, m, ...]) returns
>-           array of random numbers randomly distributed with specified mean and
>-           standard deviation"""
>-        s = standard_normal(shape)
>-        return s * std + mean
>+    """normal(mean, std, n) or normal(mean, std, [n, m, ...])
>+
>+    Returns array of random numbers randomly distributed with
>+    specified mean and standard deviation
>+    """
>+    s = standard_normal(shape)
>+    return s * std + mean
> 
> def multivariate_normal(mean, cov, shape=[]):
>-       """multivariate_normal(mean, cov) or multivariate_normal(mean, cov, [m, n, ...])
>-          returns an array containing multivariate normally distributed random numbers
>-          with specified mean and covariance.
>-
>-          mean must be a 1 dimensional array. cov must be a square two dimensional
>-          array with the same number of rows and columns as mean has elements.
>-
>-          The first form returns a single 1-D array containing a multivariate
>-          normal.
>-
>-          The second form returns an array of shape (m, n, ..., cov.getshape()[0]).
>-          In this case, output[i,j,...,:] is a 1-D array containing a multivariate
>-          normal."""
>-       # Check preconditions on arguments
>-       mean = Numeric.array(mean)
>-       cov = Numeric.array(cov)
>-       if len(mean.getshape()) != 1:
>-              raise ArgumentError, "mean must be 1 dimensional."
>-       if (len(cov.getshape()) != 2) or (cov.getshape()[0] != cov.getshape()[1]):
>-              raise ArgumentError, "cov must be 2 dimensional and square."
>-       if mean.getshape()[0] != cov.getshape()[0]:
>-              raise ArgumentError, "mean and cov must have same length."
>-       # Compute shape of output
>-       if isinstance(shape, IntType): shape = [shape]
>-       final_shape = list(shape[:])
>-       final_shape.append(mean.getshape()[0])
>-       # Create a matrix of independent standard normally distributed random
>-       # numbers. The matrix has rows with the same length as mean and as
>-       # many rows are necessary to form a matrix of shape final_shape.
>-       x = ranlib.standard_normal(Numeric.multiply.reduce(final_shape))
>-       x.setshape(Numeric.multiply.reduce(final_shape[0:len(final_shape)-1]),
>-                  mean.getshape()[0])
>-       # Transform matrix of standard normals into matrix where each row
>-       # contains multivariate normals with the desired covariance.
>-       # Compute A such that matrixmultiply(transpose(A),A) == cov.
>-       # Then the matrix products of the rows of x and A has the desired
>-       # covariance. Note that sqrt(s)*v where (u,s,v) is the singular value
>-       # decomposition of cov is such an A.
>-       (u,s,v) = LinearAlgebra2.singular_value_decomposition(cov)
>-       x = Numeric.matrixmultiply(x*Numeric.sqrt(s),v)
>-       # The rows of x now have the correct covariance but mean 0. Add
>-       # mean to each row. Then each row will have mean mean.
>-       Numeric.add(mean,x,x)
>-       x.setshape(final_shape)
>-       return x
>+    """multivariate_normal(mean, cov) or multivariate_normal(mean, cov, [m, n, ...])
>+
>+    Returns an array containing multivariate normally distributed
>+    random numbers with specified mean and covariance.
>+
>+    |mean| must be a one-dimensional array. |cov| must be a square
>+    two-dimensional array with the same number of rows and columns as
>+    |mean| has elements.
>+
>+    The first form returns a single 1-D array containing a
>+    multivariate normal.
>+
>+    The second form returns an array of shape (m, n, ...,
>+    cov.getshape()[0]). In this case, output[i,j,...,:] is a 1-D array
>+    containing a multivariate normal.
>+    """
>+    # Check preconditions on arguments
>+    mean = Numeric.array(mean)
>+    cov = Numeric.array(cov)
>+    if len(mean.getshape()) != 1:
>+        raise ArgumentError, "mean must be 1 dimensional."
>+    if (len(cov.getshape()) != 2) or (cov.getshape()[0] != cov.getshape()[1]):
>+        raise ArgumentError, "cov must be 2 dimensional and square."
>+    if mean.getshape()[0] != cov.getshape()[0]:
>+        raise ArgumentError, "mean and cov must have same length."
>+    # Compute shape of output
>+    if isinstance(shape, IntType): shape = [shape]
>+    final_shape = list(shape[:])
>+    final_shape.append(mean.getshape()[0])
>+    # Create a matrix of independent standard normally distributed random
>+    # numbers. The matrix has rows with the same length as mean and as
>+    # many rows are necessary to form a matrix of shape final_shape.
>+    x = ranlib.standard_normal(Numeric.multiply.reduce(final_shape))
>+    x.setshape(Numeric.multiply.reduce(final_shape[0:len(final_shape)-1]),
>+               mean.getshape()[0])
>+    # Transform matrix of standard normals into matrix where each row
>+    # contains multivariate normals with the desired covariance.
>+    # Compute A such that matrixmultiply(transpose(A),A) == cov.
>+    # Then the matrix products of the rows of x and A has the desired
>+    # covariance. Note that sqrt(s)*v where (u,s,v) is the singular value
>+    # decomposition of cov is such an A.
>+    (u,s,v) = LinearAlgebra2.singular_value_decomposition(cov)
>+    x = Numeric.matrixmultiply(x*Numeric.sqrt(s),v)
>+    # The rows of x now have the correct covariance but mean 0. Add
>+    # mean to each row. Then each row will have mean mean.
>+    Numeric.add(mean,x,x)
>+    x.setshape(final_shape)
>+    return x
> 
> def exponential(mean, shape=[]):
>-    """exponential(mean, n) or exponential(mean, [n, m, ...]) returns array
>-      of random numbers exponentially distributed with specified mean"""
>-   # If U is a random number uniformly distributed on [0,1], then
>-   #      -ln(U) is exponentially distributed with mean 1, and so
>-   #      -ln(U)*M is exponentially distributed with mean M.
>+    """exponential(mean, n) or exponential(mean, [n, m, ...])
>+
>+    Returns array of random numbers exponentially distributed with
>+    specified mean
>+    """
>+    # If U is a random number uniformly distributed on [0,1], then
>+    #      -ln(U) is exponentially distributed with mean 1, and so
>+    #      -ln(U)*M is exponentially distributed with mean M.
>     x = random(shape)
>     Numeric.log(x, x)
>     Numeric.subtract(0.0, x, x)
>@@ -160,52 +188,79 @@
>     return _build_random_array(ranlib.gamma, (a, r), shape)
> 
> def F(dfn, dfd, shape=[]):
>-    """F(dfn, dfd) or F(dfn, dfd, [n, m, ...]) returns array of F distributed random numbers with dfn degrees of freedom in the numerator and dfd degrees of freedom in the denominator."""
>+    """F(dfn, dfd) or F(dfn, dfd, [n, m, ...])
>+
>+    Returns array of F distributed random numbers with dfn degrees of
>+    freedom in the numerator and dfd degrees of freedom in the
>+    denominator.
>+    """
>     return ( chi_square(dfn, shape) / dfn) / ( chi_square(dfd, shape) / dfd)
> 
> def noncentral_F(dfn, dfd, nconc, shape=[]):
>-    """noncentral_F(dfn, dfd, nonc) or noncentral_F(dfn, dfd, nonc, [n, m, ...]) returns array of noncentral F distributed random numbers with dfn degrees of freedom in the numerator and dfd degrees of freedom in the denominator, and noncentrality parameter nconc."""
>+    """noncentral_F(dfn, dfd, nonc) or noncentral_F(dfn, dfd, nonc, [n, m, ...])
>+
>+    Returns array of noncentral F distributed random numbers with dfn
>+    degrees of freedom in the numerator and dfd degrees of freedom in
>+    the denominator, and noncentrality parameter nconc.
>+    """
>     return ( noncentral_chi_square(dfn, nconc, shape) / dfn ) / ( chi_square(dfd, shape) / dfd )
> 
> def chi_square(df, shape=[]):
>-    """chi_square(df) or chi_square(df, [n, m, ...]) returns array of chi squared distributed random numbers with df degrees of freedom."""
>+    """chi_square(df) or chi_square(df, [n, m, ...])
>+
>+    Returns array of chi squared distributed random numbers with df
>+    degrees of freedom.
>+    """
>     return _build_random_array(ranlib.chisquare, (df,), shape)
> 
> def noncentral_chi_square(df, nconc, shape=[]):
>-    """noncentral_chi_square(df, nconc) or chi_square(df, nconc, [n, m, ...]) returns array of noncentral chi squared distributed random numbers with df degrees of freedom and noncentrality parameter."""
>+    """noncentral_chi_square(df, nconc) or chi_square(df, nconc, [n, m, ...])
>+
>+    Returns array of noncentral chi squared distributed random numbers
>+    with df degrees of freedom and noncentrality parameter.
>+    """
>     return _build_random_array(ranlib.noncentral_chisquare, (df, nconc), shape)
> 
> def binomial(trials, p, shape=[]):
>-    """binomial(trials, p) or binomial(trials, p, [n, m, ...]) returns array of binomially distributed random integers.
>+    """binomial(trials, p) or binomial(trials, p, [n, m, ...])
>+
>+    Returns array of binomially distributed random integers.
> 
>-           trials is the number of trials in the binomial distribution.
>-           p is the probability of an event in each trial of the binomial distribution."""
>+    |trials| is the number of trials in the binomial distribution.
>+    |p| is the probability of an event in each trial of the binomial
>+    distribution.
>+    """
>     return _build_random_array(ranlib.binomial, (trials, p), shape)
> 
> def negative_binomial(trials, p, shape=[]):
>-    """negative_binomial(trials, p) or negative_binomial(trials, p, [n, m, ...]) returns
>-           array of negative binomially distributed random integers.
>+    """negative_binomial(trials, p) or negative_binomial(trials, p, [n, m, ...])
>+
>+    Returns array of negative binomially distributed random integers.
> 
>-           trials is the number of trials in the negative binomial distribution.
>-           p is the probability of an event in each trial of the negative binomial distribution."""
>+    |trials| is the number of trials in the negative binomial
>+    distribution. |p| is the probability of an event in each trial of
>+    the negative binomial distribution.
>+    """
>     return _build_random_array(ranlib.negative_binomial, (trials, p), shape)
> 
> def multinomial(trials, probs, shape=[]):
>-    """multinomial(trials, probs) or multinomial(trials, probs, [n, m, ...]) returns
>-           array of multinomial distributed integer vectors.
>+    """multinomial(trials, probs) or multinomial(trials, probs, [n, m, ...])
>+
>+    Returns array of multinomial distributed integer vectors.
>+
>+    |trials| is the number of trials in each multinomial distribution. 
>+    |probs| is a one dimensional array. There are len(prob)+1 events. 
>+    prob[i] is the probability of the i-th event, 0<=i<len(prob). The
>+    probability of event len(prob) is 1.-Numeric.sum(prob).
> 
>-           trials is the number of trials in each multinomial distribution.
>-           probs is a one dimensional array. There are len(prob)+1 events. 
>-           prob[i] is the probability of the i-th event, 0<=i<len(prob).
>-           The probability of event len(prob) is 1.-Numeric.sum(prob).
>-
>-       The first form returns a single 1-D array containing one multinomially
>-           distributed vector.
>-
>-           The second form returns an array of shape (m, n, ..., len(probs)).
>-           In this case, output[i,j,...,:] is a 1-D array containing a multinomially
>-           distributed integer 1-D array."""
>-        # Check preconditions on arguments
>+    The first form returns a single 1-D array containing one
>+    multinomially distributed vector.
>+
>+    The second form returns an array of shape (m, n, ..., len(probs)). 
>+    In this case, output[i,j,...,:] is a 1-D array containing a
>+    multinomially distributed integer 1-D array.
>+    """
>+    # Check preconditions on arguments
>     probs = Numeric.array(probs)
>     if len(probs.getshape()) != 1:
>         raise ArgumentError, "probs must be 1 dimensional."
>@@ -215,14 +270,18 @@
>     final_shape.append(probs.getshape()[0]+1)
>     x = ranlib.multinomial(trials, probs.astype(Numeric.Float32),
>                            Numeric.multiply.reduce(shape))
>-        # Change its shape to the desire one
>+    # Change its shape to the desire one
>     x.setshape(final_shape)
>     return x
> 
> def poisson(mean, shape=[]):
>-    """poisson(mean) or poisson(mean, [n, m, ...]) returns array of poisson
>-           distributed random integers with specifed mean."""
>+    """poisson(mean) or poisson(mean, [n, m, ...])
>+
>+    Returns array of poisson distributed random integers with specifed
>+    mean.
>+    """
>     return _build_random_array(ranlib.poisson, (mean,), shape)
>+
> 
> def test():
>     import test as _test
>Index: Packages/RandomArray2/Lib/__init__.py
>===================================================================
>RCS file: /cvsroot/numpy/numarray/Packages/RandomArray2/Lib/__init__.py,v
>retrieving revision 1.1
>diff -u -r1.1 __init__.py
>--- Packages/RandomArray2/Lib/__init__.py	21 Jun 2002 18:25:29 -0000	1.1
>+++ Packages/RandomArray2/Lib/__init__.py	22 Aug 2002 22:11:16 -0000
>@@ -2,3 +2,6 @@
> 
> from RandomArray2 import *
> 
>+__doc__ = RandomArray2.__doc__ + """
>+See RandomArray2.RandomArray2 for more information.
>+"""
>


-- 
Todd Miller 			jmiller at stsci.edu
STSCI / SSG