[SciPy-user] Simple combinatorics with Numpy

T J tjhnson@gmail....
Tue Oct 21 16:25:49 CDT 2008

On Tue, Oct 21, 2008 at 2:52 AM, M F wrote:
> Imagine I have an array x=([1,3,5,8]) from which I draw, with replacement, a
> sample of, say, 3 numbers. The ordering of the sample is unimportant. I want to
> calculate the histogram of the averages of the samples (i.e., avg(1,1,1)=1,
> avg(1,3,5)=3, etc). For that, I need to *enumerate* all the possible unordered
> samples of size 3 that can be drawn from x, i.e.,

There was a post in numpy-discussion a couple of days ago which
provided a possible solution:


def boxings(n, k):
    seq, i = [n]*k + [0], k
    while i:
        yield tuple(seq[i] - seq[i+1] for i in xrange(k))
        i = seq.index(0) - 1
        seq[i:k] = [seq[i] - 1] * (k-i)

This enumerates the ways to place n indistinguishable tokens into k
distinguishable boxes (order of boxes is not of concern) when there is
no maximum on the occupancy of any box.  You want to make j unordered
selections from m distinguishable items when the items are replaced
after each selection.

n <-> j
k <-> m

In the example you suggested, boxings(3,4) does the trick, but this
gives you the number of times each item was selected.  There is
probably a simple way to directly generate the samples, but it is
possible to do it using the above function too.

def samples_ur(items, k):
    """Returns k unordered samples (with replacement) from items."""
    n = len(items)
    for sample in boxings(k, n):
        selections = [[items[i]]*count for i,count in enumerate(sample)]
        yield tuple([x for sel in selections for x in sel])

>From there, you should be able to do what you want (eg average).

> and then, I need to compute, for each of them, the weight c/(4**3), where c is a
> multinomial factor that takes into account possible repetitions of elements in
> the sample.

What exactly are you looking for with c?  Do you want to know the
number of permutations for each enumerated item?  This isn't the
prettiest code, but if so, the attached file does what you want.
-------------- next part --------------
A non-text attachment was scrubbed...
Name: samples.py
Type: text/x-python
Size: 1139 bytes
Desc: not available
Url : http://projects.scipy.org/pipermail/scipy-user/attachments/20081021/ec29db9d/attachment.py 

More information about the SciPy-user mailing list