Random permutations
Mok-Kong Shen
Hi,
If I don't err, one generally employs Algorithm P in Knuth's book
(vol.2 3rd ed. p.145) to produce random permutations of t objects
X_1, X_2, .... X_t. Would the algorithm below also be acceptable?
Do t times the following:
Get 2 integer random numbers i, and j, uniformly distributed
in [1, t]. Exchange X_i and X_j.
Thanks,
M. K. Shen
Richard Heathfield
Model it, and see what kind of distribution you get.
(I would go into a lot more detail, but your name is familiar to me,
and I *know* you're bright enough to know exactly what I mean and
exactly how to carry it out.)
--
Richard Heathfield <http://www.cpax.org.uk>
Email: -http://www. +rjh@
"Usenet is a strange place" - dmr 29 July 1999
Sig line vacant - apply within
Pascal J. Bourguignon
A precise statistical analysis would be required, but AFAICUI this
doesn't result in the same distribution of resulting permutations,
since you can get two or more times the same set {i,j}, and therefore
cancel some exchanges.
--
__Pascal Bourguignon__
Mok-Kong Shen
> Mok-Kong Shen wrote:
> ...... produce random permutations of t objects
>> X_1, X_2, .... X_t. Would the algorithm below also be acceptable?
>>
>> Do t times the following:
>>
>> Get 2 integer random numbers i, and j, uniformly distributed
>> in [1, t]. Exchange X_i and X_j.
>
> Model it, and see what kind of distribution you get.
>
> (I would go into a lot more detail, but your name is familiar to me,
> and I *know* you're bright enough to know exactly what I mean and
> exactly how to carry it out.)
Thank you but you overestimated my capability. I am not a mathematician
and have really only superficial and meager knowledge in probability
and statistics. I sould be very grateful, if you would kindly say
something more on how to do the modelling or point to 'concrete'
directions (with literatuere) that I could first attempt to follow
myself and eventually later ask further help from you, in case I
don't fare very far.
BTW, my reasoning was as follows: Each permutation is the result of
a sequence of exchanges of 2 elements. Doing a number of random
exchanges doesn't bias any particular element and should thus result
in a random permutation.
M. K. Shen
Pascal J. Bourguignon
Let's take an example with 3 elements.
An equiprobable distribution of the permutations would be:
abc 1/6
acb 1/6
bac 1/6
bca 1/6
cab 1/6
cba 1/6
Now, applying your algorithm (0-based as it should) we get:
1/9 1/9 1/9 1/9 1/9 1/9 1/9 1/9 1/9
(0 0) (0 1) (0 2) (1 0) (1 1) (1 2) (2 0) (2 1) (2 2)
abc --> abc bac cba bac abc acb cba acb abc
or: abc 3/9
acb 2/9
bac 2/9
bca 0/9
cab 0/9
cba 2/9
Applying it again at step 2:
abc 1/3
abc 1/3
acb 2/9
bac 2/9
bca 0
cab 0
cba 2/9
acb 2/9
acb 1/3
abc 2/9
cab 2/9
cba 0
bac 0
bca 2/9
bac 2/9
bac 1/3
bca 2/9
abc 2/9
acb 0
cba 0
cab 2/9
bca 0/9
cab 0/9
cba 2/9
cba 1/3
cab 2/9
bca 2/9
bac 0
acb 0
abc 2/9
or: abc 1/3*1/3 + 2/9*2/9 + 2/9*2/9 + 0 + 0 + 2/9*2/9 = 7/27
acb 1/3*2/9 + 2/9*1/3 + 0 + 0 + 0 + 0 = 4/27
bac 1/3*2/9 + 0 + 2/9*1/3 + 0 + 0 + 0 = 4/27
bca 0 + 2/9*2/9 + 2/9*2/9 + 0 + 0 + 2/9*2/9 = 4/27
cab 0 + 2/9*2/9 + 2/9*2/9 + 0 + 0 + 2/9*2/9 = 4/27
cba 1/3*2/9 + 0 + 0 + 0 + 0 + 2/9*1/3 = 4/27
Applying it a third time:
abc 7/27
abc 1/3
acb 2/9
bac 2/9
bca 0
cab 0
cba 2/9
acb 4/27
acb 1/3
abc 2/9
cab 2/9
cba 0
bac 0
bca 2/9
bac 4/27
bac 1/3
bca 2/9
abc 2/9
acb 0
cba 0
cab 2/9
bca 4/27
bca 1/3
bac 2/9
cba 2/9
cab 0
abc 0
acb 2/9
cab 4/27
cab 1/3
cba 2/9
acb 2/9
abc 0
bca 0
bac 2/9
cba 4/27
cba 1/3
cab 2/9
bca 2/9
bac 0
acb 0
abc 2/9
which gives:
abc 7/27*1/3 + 4/27*2/9 + 4/27*2/9 + 0 + 0 + 4/27*2/9 = 5/27
acb 7/27*2/9 + 4/27*1/3 + 0 + 4/27*2/9 + 4/27*2/9 + 0 = 14/81
bac 7/27*2/9 + 0 + 4/27*1/3 + 4/27*2/9 + 4/27*2/9 + 0 = 14/81
bca 0 + 4/27*2/9 + 4/27*2/9 + 4/27*1/3 + 0 + 4/27*2/9 = 4/27
cab 0 + 4/27*2/9 + 4/27*2/9 + 0 + 4/27*1/3 + 4/27*2/9 = 4/27
cba 7/27*2/9 + 0 + 0 + 4/27*2/9 + 4/27*2/9 + 4/27*1/9 = 14/81
Clearly this is not an equiprobable distribution of the permutations.
--
__Pascal Bourguignon__
James Dow Allen
> Would the algorithm below [produce random t-sized permutations]?
>
> Do t times the following:
> Get 2 integer random numbers i, and j, uniformly distributed
> in [1, t]. Exchange X_i and X_j.
No. One way to see this is might be to do it k times instead of t
times. If k is very small (say k=1 or even k=0) it obviously
won't work. As k increases, the randomness gets better *but* there's
nothing special about k=t that would suddenly make it perfectly
random. (A Principle of Continuity applies.) In this method perfect
(random) uniformity would be approached only as k approaches infinity.
Perhaps another way to see this is to note you're picking only
2t numbers from the set of size t, so overlooking a number entirely
won't be uncommon, and that prevents the randomness you seek.
BTW, the Algo P you refer to involves selecting t-1 random integers
from 1...t, 1...t-1, etc. If you're worried about the time
wasting of multiple random() invocations, you can get several
numbers from one call by doing appropriate /'s and %'s.
Somone whose posts usually have much wisdom wrote:
> Model it, and see what kind of distribution you get.
IMO, proofs (or intuitions) of mathematical fact, when available,
are superior to experimentation.
James
Richard Heathfield
<40fb2ac7-4e11-4e5d...@m33g2000pri.googlegroups.com>,
James Dow Allen wrote:
<snip>
> Somone whose posts usually have much wisdom wrote:
>> Model it, and see what kind of distribution you get.
s/usually/always/ :-)
> IMO, proofs (or intuitions) of mathematical fact, when available,
> are superior to experimentation.
Sometimes, however, the fastest way to disprove a hypothesis is to
acquire a counter-example, which you get by experimentation. In the
case of FLT, that wouldn't be bright! But in this case, even a simple
model will quickly throw up a counter-example, as Pascal B
demonstrated (which is why I didn't bother). I note that his
counter-example appeared somewhat earlier than your mathematical
insights, supporting my hypothesis that experimentation can sometimes
be faster than a more analytical approach.
That is most emphatically *not* to say, however, that an analytical
approach has no value, or that experimentation is always the best
option.
James Dow Allen
> I note that his
> counter-example appeared somewhat earlier than your mathematical
> insights, supporting my hypothesis that experimentation can sometimes
> be faster than a more analytical approach.
It almost sounds like you think I'm glued to my Internet
connection (like some of you :-), but it took me some hours
to come up with this mathematical insight!
If it does seem like I'm posting more often these last few
days, it's because for the first time in my life ... I have
Internet at home!!! (My daughter wanted it for her birthday.
My prediction -- that it would aggravate her father's insomnia --
is proving correct.)
The connection is slow (115kb/s max) and bugs appear. One
I hadn't seen for a while has reappeared: Google Groups
will not let me see the ends of Messages 10, 20, 30, etc.
in "Threaded mode"! I'll look for a Google-hack to get around
that, but meanwhile please don't post *good* messages at
those ordinal positions. :-)
James
pete
I use this algorithm:
#define LU_RAND(S) ((S) * 69069 + 362437)
long unsigned shuffle(long unsigned *array,
long unsigned n,
long unsigned seed)
{
long unsigned i, r;
array[0] = 0;
for (i = 1; n > i; ++i) {
seed = LU_RAND(seed);
r = seed % (i + 1);
array[i] = 0;
array[i] = array[r];
array[r] = i;
}
return seed;
}
--
pete
Patricia Shanahan
Is it just curiosity, or do you have some reason to avoid using the
normal, provably correct algorithm? If you do need to avoid the normal
algorithm, there may be some other way of dealing with the issue.
Patricia
mike
online.de says...
For any randomly permutated list of t elements, the probability that any
chosen element (for example the original first element X_1) is in its
original position is precisely 1/t.
For your algorithm, if an item is not chosen in any of the t swaps, it
will remain in its first position. So the probability it remains unmoved
is at least equal to ((t-1)/t)^2t. For larger t, this limits to 1/e.
Note I said that the probability is at least this, as there is also the
probability that an element is switched out and then switched in again.
So clearly the statistical distribution of the swap permutation will
differ from that of the random permutation as 1/t .NE. 1/e
You could carry out more swaps and get closer to a 'random' distribution
but
a) it would be a very inefficient process as you would need at least
t*ln(t) swaps before you began to approach a random distribution, and
b) it would still not be exactly random.
Mike
Dave Searles
> In
> <40fb2ac7-4e11-4e5d...@m33g2000pri.googlegroups.com>,
> James Dow Allen wrote:
>
> <snip>
>
>> Somone whose posts usually have much wisdom wrote:
>>> Model it, and see what kind of distribution you get.
>
> s/usually/always/ :-)
>
>> IMO, proofs (or intuitions) of mathematical fact, when available,
>> are superior to experimentation.
>
> Sometimes, however, the fastest way to disprove a hypothesis is to
> acquire a counter-example, which you get by experimentation. In the
> case of FLT, that wouldn't be bright! But in this case, even a simple
> model will quickly throw up a counter-example, as Pascal B
> demonstrated (which is why I didn't bother). I note that his
> counter-example appeared somewhat earlier than your mathematical
> insights, supporting my hypothesis that experimentation can sometimes
> be faster than a more analytical approach.
>
> That is most emphatically *not* to say, however, that an analytical
> approach has no value, or that experimentation is always the best
> option.
In fact, it seems that an analytical approach has value in proving
something true, and experimentation has value in proving something false.
Mok-Kong Shen
[snip]
> BTW, the Algo P you refer to involves selecting t-1 random integers
> from 1...t, 1...t-1, etc. If you're worried about the time
> wasting of multiple random() invocations, you can get several
> numbers from one call by doing appropriate /'s and %'s.
I have learnt from several follow-ups that my scheme doesn't
work correctly. On the other hand, I don't understand your last
sentence above. Anyway, I doubt the usefulness of % here, since,
if t is uniformly distributed in [0, n-1], then for m<n the
values t mod m wouldn't be uniformly distributed in [0, m-1]
in general, if I don't err.
M. K. Shen
Pascal J. Bourguignon
You err, in the same way as your algorithm was wrong.
Let's take an example, 11 % 3:
. . .
. . .
. . .
. .
----------------
4/11 4/11 3/11
t mod m can be uniformly distributed only if n is a multiple of m.
--
__Pascal Bourguignon__
Mok-Kong Shen
I suppose you misread what I wrote. I wrote: "... t mod m WOULDN'T
be uniformly distributed in general.
M. K. Shen
James Dow Allen
> James Dow Allen wrote:
>
> [snip]
>
> > BTW, the Algo P you refer to involves selecting t-1 random integers
> > from 1...t, 1...t-1, etc. If you're worried about the time
> > wasting of multiple random() invocations, you can get several
> > numbers from one call by doing appropriate /'s and %'s.
>
Perhaps I should have elaborated. Suppose you want three
random numbers, 0 <= a < 10, 0 <= b < 9, 0 <= c < 8.
Call random() once to get a random number 0 <= z < 10*9*8.
Then do
a = z % 10;
b = (z / 10) % 9;
c = z / 90;
Obviously you can get more than 3 random numbers from a
single random() this way, as long as the product of their
ranges is compatible with your ranged_integer_from_random()
mechanism. I *won't* give a ranged_integer_from_random()
implementation here, since that topic has been beaten
to death over and over.
James Dow Allen
Pascal J. Bourguignon
Mok-Kong Shen
But a follow-up of Bourguignon shows that using % wouldn't
give a 'really' random number in general. See his post.
M. K. Shen
James Dow Allen
> James Dow Allen schrieb:
> > Perhaps I should have elaborated. Suppose you want three
> > random numbers, 0 <= a < 10, 0 <= b < 9, 0 <= c < 8.
> > Call random() once to get a random number 0 <= z < 10*9*8.
> > Then do
> > a = z % 10;
> > b = (z / 10) % 9;
> > c = z / 90;
> > Obviously you can get more than 3 random numbers from a
> > single random() this way, as long as the product of their
> > ranges is compatible with your ranged_integer_from_random()
> > mechanism. I *won't* give a ranged_integer_from_random()
> > implementation here, since that topic has been beaten
> > to death over and over.
>
> But a follow-up of Bourguignon shows that using % wouldn't
> give a 'really' random number in general. See his post.
Wrong again.
Pascal wrote:
> if t is uniformly distributed in [0, n-1], then ...
> t mod m can be uniformly distributed only if n is a multiple of m.
I wrote:
> > Call random() once to get a random number 0 <= z < 10*9*8.
> > a = z % 10;
"n" is 10*9*8. n is a multiple of 10.
> > b = (z / 10) % 9;
Since "n" was a multiple of 10, (z / 10) is uniform on 1 ... 9*8-1
That is, the *new* "n" is 9*8
> > b = (z / 10) % 9;
"n" is 9*8. n is a multiple of 9.
> > c = z / 90;
Since z is a uniform variate on 0 ... 10*9*8-1
"n" is now 10*9*8 as before. n is a multiple of 90.
Are we getting clear yet?
James Dow Allen
Sorry for being brusque! (I don't mind doing someone's
homework for free and being wrongly contradicted about it,
but this isn't a chatroom, and the person whose homework
is done for free is expected to at least try to put his
own thinking cap on before repeating the same wrong
contradiction three times!)
On Oct 21, 3:17 am, James Dow Allen <jdallen2...@yahoo.com> wrote:
> Suppose you want three
> random numbers, 0 <= a < 10, 0 <= b < 9, 0 <= c < 8.
> Call random() once to get a random number 0 <= z < 10*9*8.
The above itself may be non-trivial. But you'd have this
problem (developing a uniform variate whose range is not a
power of two) with or without the technique which is the
subject of this subthread.
As Pascal's little diagram shows, you want "n" which is a
multiple of your range(s). The 10*9*8 works out just right!
> Then do
> a = z % 10;
> b = (z / 10) % 9;
> c = z / 90;
I don't want to be a dictator! Let's vote!
The claim is that
IF z is uniform with range 8*9*10
THEN a, b, c (as developed in the above code fragment)
are uniform (and independent of each other)
with ranges 10, 9, 8 respectively.
I vote YES. Pascal, are you awake?
(I'd be *very* surprised if this technique isn't *very*
familiar to comp.programmer's who work with such
randomizations ... but then I've been surprised before.)
James
Pascal J. Bourguignon
This is correct. But as you said, it all depends on a
ranged_integer_from_random() function. Using z=random() would be bad.
And once you have a ranged_integer_from_random(), you don't need to
call it on p*q*r, since calling it on p, then on q, then on r will
essentially do the same thing you do here with / and %.
--
__Pascal Bourguignon__
Pascal J. Bourguignon
The claim is correct.
What is not correct was the "call random() once" premise.
> (I'd be *very* surprised if this technique isn't *very*
> familiar to comp.programmer's who work with such
> randomizations ... but then I've been surprised before.)
>
> James
--
__Pascal Bourguignon__
Mok-Kong Shen
Several people in this thread have kindly shown (and with different
arguments) that the scheme quoted above cannot work correctly.
There was indeed, as you questioned, a motivation that (somehow
indirectly) led to my original post. The scenario I consider can be
stated with the following example:
Given a random (8 bit) byte sequence U_i, one desires to use it to
generate random permutations of 256 objects numbered 0..255. In
applying Algorithm P in Knuth's book, one needs at each step a random
number in [0,m], where m is a value decreasing from 255 by 1 at each
processing step. Now if one discards U_i as long as U_i > m, i.e.
waiting for the first U_i to appear that is less or equal to m, one
would evidently waste in general quite a large amount of U_i, when
m becomes small. This is certainly very undesirable. If on the other
hand one uses R = U_i mod m, then R isn't uniformly distributed in
[0,m].
I suppose that the best one could do in the given situation is
apparently the following (rather trivial) solution:
For an m above, let k be such that 2^k > m >= 2^(k-1), Take from
U_i k bits and save the remaining 8-k bits for future use (i.e. to
be concatenated with the bits of the following element of the byte
sequence). We have now a k bit random variable V_i > m such that the
probability of V_i <= m is at least 0.5. That is, we have in total
a waste of the given random material less than 50 percent.
M. K. Shen
Mok-Kong Shen
[skip]
> sequence). We have now a k bit random variable V_i > m such that the
> probability of V_i <= m is at least 0.5. .........
Sorry, typo. Please read:
.......... We have now a k bit random variable V_i such that the
probability of V_i <= m is at least 0.5. .........
M. K. Shen
James Dow Allen
wrote:
> James Dow Allen <jdallen2...@yahoo.com> writes:
> > The claim is that
> > IF z is uniform with range 8*9*10
> > are uniform (and independent of each other)
> > with ranges 10, 9, 8 respectively.
>
> > I vote YES. Pascal, are you awake?
>
> The claim is correct.
Well, thanks! The vote is now 2-0!
> What is not correct was the "call random() once" premise.
That's right, as I'd already implied, albeit tersely.
I *should* have written "approximately once"
instead of "once"! :-)
But *that* issue was not germane to this particular
subthread and has certainly been discussed here recently.
James
James Dow Allen
> There was indeed, as you questioned, a motivation that (somehow
> indirectly) led to my original post. The scenario I consider can be
> stated with the following example:
> generate random permutations of 256 objects numbered 0..255.
> would evidently waste in general quite a large amount of U_i, when
> m becomes small.
So *this* was your problem the whole time! We *might* have so
guessed had the "t" in your first post been an explicit power of two.
As I just wrote upthread the problem and its solution is frequently
discussed here. Normal solution is to acquire a 31-bit random
integer x; take x's quotient q AND remainder r, both modulo t.
IF (q+1) * t <= 2^31, THEN r is your desired uniform t-range
variate; ELSE acquire a new 31-bit random integer and try again.
Hope this helps.
You are correct to worry about "wasted" random bits, at least if
you're perfectionistic. There are *two* sources of such waste here.
(1) If t is much smaller than 2^31 you're wasting bits. For example,
if t = 2^8, you're wasting 23 bits each time (assuming 31-bit
generator).
I showed you already an easy way to avoid this.
(2) When you "acquire a new random and retry" in the procedure above,
the rejected random number is partially wasted. The loss, however,
will be trivially small unless t is slightly larger than a large
power of two *and* since, almost by definition, most of the high-order
bits on such a rejected number are all ones, recovering the "still
truly
random" but otherwise wasted bits is non-trivial.
Hope this helps too.
James
James Dow Allen
wrote:
> >> > mechanism. I *won't* give a ranged_integer_from_random()
> >> > implementation here, since that topic has been beaten
> >> > to death over and over.
>
> call it on p*q*r, since calling it on p, then on q, then on r will
> essentially do the same thing you do here with / and %.
I just saw this message. Gosh. As I explained a few moments ago,
there are at least two (2) issues. Since no one has written the
man page for ranged_integer_from_random() yet, we don't know
whether it saves otherwise wasted bits (as your message implies that
you think it should) or is a small backend that can be used for
such a random-bit-saving funxtion. I was showing the OP the
simple way to do that bit-saving (or offering a glimpse at how to
*write* the "better" ranged_integer_from_random(), if you wish).
These threads would become less exasperating if knowledgeable
posters who *do not disagree* with each other would refrain from
misconstructions whose sole purpose is to find a way to quibble.
James
Mok-Kong Shen
[snip]
> .... If on the other hand one uses R = U_i mod m, then R isn't
Pardon, again typo. Please read R = U_i mod (m+1).
M. K. Shen
Mok-Kong Shen
> (1) If t is much smaller than 2^31 you're wasting bits. For example,
> if t = 2^8, you're wasting 23 bits each time (assuming 31-bit
> generator).
> I showed you already an easy way to avoid this.
A not uncommon situation is that one has a good quality random bit
sequence (in byte or word format) at hand and it seems preferable
(or even desired) to use that directly than using any function to
compute something like random().
> (2) When you "acquire a new random and retry" in the procedure above,
> the rejected random number is partially wasted. The loss, however,
> will be trivially small unless t is slightly larger than a large
> power of two *and* since, almost by definition, most of the high-order
> bits on such a rejected number are all ones, recovering the "still
> truly
> random" but otherwise wasted bits is non-trivial.
With Knuth's Algoithm P, one would often have the situation of a
decreasing t that approaches a power of 2. That problem couldn't
naturally be avoided, unfortunately.
Thanks,
M. K. Shen
James Dow Allen
> James Dow Allen wrote:
> > ... The loss, however,
> > will be trivially small unless t is slightly larger than a large
> With Knuth's Algoithm P, one would often have the situation of a
> decreasing t that approaches a power of 2. That problem couldn't
> naturally be avoided, unfortunately.
I wrote:
LARGE power of two; LARGE power of two; LARGE power of two!
Let me repeat that:
LARGE power of two; LARGE power of two; LARGE power of two!
LARGE power of two; LARGE power of two; LARGE power of two!
LARGE power of two; LARGE power of two; LARGE power of two!
LARGE power of two; LARGE power of two; LARGE power of two!
LARGE power of two; LARGE power of two; LARGE power of two!
For a number that is slightly larger than a SMALL power of two,
for example 2055, the procedure I described would require
you to loop less than one time in a million trials, on average.
If you ARE shuffling a TRULY immense deck, perhaps you need
a high-level redesign anyway.
Hope this helps. If not, I give up.
James
Mok-Kong Shen
> I wrote:
> LARGE power of two; LARGE power of two; LARGE power of two!
I am afraid that I yet misunderstand you. What I mean is this: I am
interested in permutations of number of objects up to, say, 1024.
Now, with Algorithm P, one starts using t of 10 bits for a random
value m, whose max size decreases in each processing step by 1 till
the max size bedomes 511, where one begins using t of 9 bits etc.
Does that conform to your idea of applying Algorithm P or have
I totally misunderstood you? (Eventually my apology in advance.)
Thanks,
M. K. Shen
James Dow Allen
> I am afraid that I yet misunderstand you...
I will answer your question, but first I have one of my own.
Throughout this "conversation" have you been thinking that
James is confused or that Mok-Kong is confused? (or both? :-)
I ask this because IF you assumed I knew what I was talking
about AND cleared away your preconceptions to just imagine
the implementation I describe, THEN I believe an intelligent
programmer (which I assume you are) would have been able to
figure this out yesterday. Instead you're working with
preconceptions different from mine, and listening only to
bits and pieces of my story, trying to fit them to your
preconception. Read my answer carefully; write actual
working code if necessary until you grasp the principles.
Meanwhile please answer the question above.
> Now, with Algorithm P, one starts using t of 10 bits for a random
> value m,...
Now, let's next agree on what YOUR question really is.
The procedure you just described to shuffle a 1024-sized
set WORKS FINE, but suffers from the following defect:
When you need a variate with range 513 in your procedure,
you take "10 random bits" (that is, a random number
0 <= x < 1024). Then, to ensure that the uniform
distribution you seek is emulated EXACTLY, you DISCARD
that random number whenever x >= 513 and get ANOTHER
random number.
In that way, your procedure is CORRECT but spends on average
19.98 random bits to select the 513-range variate, when
a mere 9.02 bits would be sufficient in theory. You have
some reason (beyond the scope of our conversation) for
avoiding this wastage of random bits. You want to know how
to approach a theoretical optimum closer to 9.02 bits,
rather than the wastage implicit in the "loop until x < 513"
mechanism. (In the sequel I write "x < 513" as shorthand
for "x < whatever"; I hope *that* doesn't cause confusion.)
Is this much correct? Do I understand your question?
Assuming I *do* understand your question; let me first
say:
Since a typical random() function returns 31
bits, I'll guess that 99% of programmers doing a shuffle
use 31 bits for EACH variate. This means that in the
specific case we consider, they WASTE almost 22 bits!
Congratulations on being perfectionistic enough to want
to avoid this; BUT, if your special need is to avoid that
waste, don't you think it would have clarified this
discussion to MENTION it? (Or am I wrong, and the problem
is that you just didn't think of the "loop until x < 513"
way of perfecting the distribution???)
Now I'll assume you KNOW how to shuffle correctly and
are ONLY concerned about how to minimize the wastage
of your random bits!!
There are two obvious approaches to do this.
Approach A:
A is based on what I showed you earlier. I'll work a
specific example and leave it as an exercise to code
a general routine. For the specific example, we need
variates of ranges 11, 10 and 9. These are each slightly
larger than 2^3 so not only do we need 4+4+4 bits
minimum to develop them in the naive way, we run
a high "wastage risk" and will spend an average of
about 21 random bits altogether for the three variates.
HOWEVER, because 11*10*9 = 990 is slightly LESS than
the power of two 2^10, you can develop the variates
using the technique based on "/" and "%" from only
10 bits (plus wastage of about 0.1 bits).
Does this make sense? Write and test a program
right now to develop these three variates just as
described. To fully use this technique you'd
want GENERAL code rather than just the 11,10,9
example but we'll leave that as a future exercise....
Approach B uses a different principle altogether.
You do things the old naive way BUT whenever you
REJECT a random number (because of the inevitable
"loop until x < 513" mechanism) you SAVE some of
that otherwise-discarded random number! BECAUSE
the portion of the random number that can be SAVED
will NOT be a whole number of bits, this approach
is NON-TRIVIAL. I seem to be in a patient mood,
but I WON'T show you how to do this because
(a) Approach A works already.
(b) I think it's a better-than-even chance
I've lost you already.
REALLY hope this helps.
James Dow Allen
Mok-Kong Shen
[snip]
I readily concede that I am layman with very poor knowledge (and IQ).
But let me stress that my assumption is that there are readily
availble very good random bits (though it is desired to use them
economically). Now, e.g. for a random number in the range [0,255],
why should one consider to write a program to do that instead of
simply taking 8 bits right away from the given bit sequence? This
explains perhaps that I am too lazy yet to consider any solution
that needs programming work. BTW, note that, since the random bits
are assumed to be of high quality, I rather doubt that any random
number generated with some good programming work would lead to much
superior results in comparison.
A counter argument could be: Where does the entropy of random()
come from? Note that the assumption is that there are very good
random bits availble, i.e. with almost full entropy in them.
For e.g. each 8 bits IS already a random number in [0,255],
without needing any programming work. If random() could continue
to deliver bits of full entropy without further ado and without
limit, that would certainly be absolutely excellent. But I don't
think this can be so, before you convince me otherwise.
About the issue of programming, please see my point above.
Thanks,
M. K. Shen