On Tue, Mar 29, 2011 at 1:29 PM, eat wrote:

Hi,

On Tue, Mar 29, 2011 at 12:00 PM, Alex Ter-Sarkissov wrote:

...

If I want to generate a string of random bits with equal probability I run

random.randint(0,2,size).

What if I want a specific proportion of bits? In other words, each bit is 1 with probability p<1/2 and 0 with probability q=1-p?

Would it be sufficient to: In []: bs= ones(1e6, dtype= int) In []: bs[randint(0, 1e6, 1e5)]= 0 In []: bs.sum()/ 1e6 Out[]: 0.904706

Or: In []: bs= ones(1e6, dtype= int) In []: bs[rand(1e6)< 1./ 10]= 0 In []: bs.sum()/ 1e6 Out[]: 0.89983

Regards, eat

...

thanks

_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

Den 29.03.2011 14:56, skrev Sturla Molden:

import numpy as np def randombits(n, p): b = (np.random.rand(n*8).reshape((n,8))< p).astype(int) return (b<< range(8)).sum(axis=1).astype(np.uint8)

n is the number of bytes. If you prefer to count in bits: def randombits(n, p): assert(n%8 == 0) b = (np.random.rand(n).reshape((n>>3,8))< p).astype(int) return (b<< range(8)).sum(axis=1).astype(np.uint8) Sturla

Den 29.03.2011 15:46, skrev Daniel Lepage:

x = (np.random.random(size)< p)

This will not work. A boolean array is not compactly stored, but an array of bytes. Only the first bit 0 is 1 with probability p, bits 1 to 7 bits are 1 with probability 0. We thus have to do this 8 times for each byte, shift left by range(8), and combine them with binary or. Also the main use of random bits is crypto, which requires the use of /dev/urandom or CrypGenRandom instead of Mersenne Twister (np.random.rand). Sturla

On Tue, Mar 29, 2011 at 09:49, Sturla Molden wrote:

Den 29.03.2011 15:46, skrev Daniel Lepage:

...

x = (np.random.random(size)<  p)

This will not work. A boolean array is not compactly stored, but an array of bytes. Only the first bit 0 is 1 with probability p, bits 1 to 7 bits are 1 with probability 0. We thus have to do this 8 times for each byte, shift left by range(8), and combine them with binary or.

It's not clear that the OP really meant bits rather than just bools. Judging by the example that he tried first, it's likely that he just wants bools (or even just 0s and 1s) and not a real string of bits compacted into bytes. -- Robert Kern "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

On Tue, Mar 29, 2011 at 11:59 AM, Sturla Molden wrote:

Den 29.03.2011 16:49, skrev Sturla Molden:

...

This will not work. A boolean array is not compactly stored, but an array of bytes. Only the first bit 0 is 1 with probability p, bits 1 to 7 bits are 1 with probability 0. We thus have to do this 8 times for each byte, shift left by range(8), and combine them with binary or. Also the main use of random bits is crypto, which requires the use of /dev/urandom or CrypGenRandom instead of Mersenne Twister

(np.random.rand).

Here's a cleaned-up one for those who might be interested :-)

How about adding it to http://www.scipy.org/Cookbook? Warren

Sturla

import numpy as np import os

def randombits(n, p, intention='numerical'): """ Returns an array with packed bits drawn from n Bernoulli trials with successrate p. """ assert (intention in ('numerical','crypto')) # number of bytes m = (n >> 3) + 1 if n % 8 else n >> 3 if intention == 'numerical': # Mersenne Twister rflt = np.random.rand(m*8) else: # /dev/urandom on Linux, Apple, et al., # CryptGenRandom on Windows rflt = np.frombuffer(os.urandom(m*64),dtype=np.uint64) rflt = rflt / float(2**64) b = (rflt.reshape((m,8))> (m*8 - n)) return b

def bitcount(a, bytewise=False): """ Count the number of set bits in an array of np.uint8. """ assert(a.dtype == np.uint8) c = a[:,np.newaxis].repeat(8, axis=1) >> range(8) & 0x01 return (c.sum(axis=1) if bytewise else c.sum())

if __name__ == '__main__':

b = randombits(int(1e6), .1, intent='numerical') print bitcount(b) # should be close to 100000 b = randombits(int(1e6), .1, intent='crypto') print bitcount(b) # should be close to 100000

_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion