When does randomness enter a shuffled "deck" of coin tosses.
2 views (last 30 days)
Show older comments
george simpson
on 12 Jul 2020
Commented: george simpson
on 17 Jul 2020
%shuffle a clean deck of 1e6 tails and heads
p=ones(1,500000); m=-ones(1,500000);C=[p ; m];a=C(:)';
%1st shuffle
A=a(1:500000);B=a(500001:1e6);C=[A ; B];C1=C(:)';
% next 1e3 shuffles
A=C1(1:5e5);B=C1(5e5+1:1e6);C=[A;B];C2=C(:)';
A=C2(1:5e5);B=C2(5e5+1:1e6);C=[A;B];C3=C(:)';
A=C3(1:5e5);B=C3(5e5+1:1e6);C=[A;B];C4=C(:)';
...
Brainlock :
How to use output of a command as input in following repeated commands
3 Comments
John D'Errico
on 12 Jul 2020
Edited: John D'Errico
on 12 Jul 2020
Almost impossible to understnd your question. My guess is you are starting with an unshuffled binary deck, of length 1e6. So half are zero, half are 1. Now you want to do a perfect riffle shuffle? (EXPLAIN WHAT YOU ARE DOING! It is difficult to help you if we must guess.)
It looks like you are effectively just cutting the deck, midway, then combineing the two halves into two rows of an array, and the stringing it out again into a row array. So in fact, that will perform a true riffle shuffle for a vector of even length. In fact, that is a good and efficient way to perform a riffle shuffle.
So, is your goal to do multiple riffle shuffles, then eventually looking to see when the result "appears" random?
A problem is it is difficult to know if a sequence is truly random. But exactly what are you trying to do? What is the goal? Of course, you do need to know that the sequence you get will always be perfectly repeatable, so never truly random. You will get always the same sequence after n such shuffles.
Accepted Answer
John D'Errico
on 12 Jul 2020
Edited: John D'Errico
on 12 Jul 2020
I put together this simple script. In it you will find perhaps the most efficient scheme for a riffle shuffle that I know of, since it can be done as a simple element permutation, just a re-indexing of the elements.
With N = 5e5, thus 1e6 total elements, you should see how fast it is.
timeit(@() B(riffleperm))
ans =
0.005221
So it would appear I can do nearly 200 such shuffles in one second.
Anyway, read down to the end, and then think about what you see.
%% Will a sequence of riffle shuffles ever produce a truly random sequence? Well, yes. And no.
% Consider a sequence of random bits, N zeros, followed by N ones. Now
% perform a perfect riffle shuffle. Repeat as often as you wish. The final
% sequence of bits after S such perfect riffle shuffles is perfectly
% deterministic, so not at all random. It may appear to be pseudo-random.
% It may approach the appearance of randomness at some level. We can test
% for things like that.
N = 5e5;
B0 = floor((0:2*N-1)/N);
% riffleperm has the property that B(riffleperm) will be a perfect riffle
% shuffle.
riffleperm = zeros(1,2*N);
riffleperm(1:2:end) = 1:N;
riffleperm(2:2:end) = N+1:2*N;
%% Now, perform S such riffleshuffles. For example, with S = 10, we might do this:
S = 10;
B = B0;
for ns = 1:S
B = B(riffleperm);
end
%% How random is B, after 10 such shuffles?
% One simple test might be to ask how many 1s are in the first half of the
% sequence. While a perfectly random sequence might seem to have exactly
% 1/2 of the 1's that would be wrong. In fact, the distribution of the
% number of 1's in either half of a truly random sequence would follow a
% binomial distribution.
sum(B(1:N))
ans =
250000
% In fact, that is a highly unlikely event. It is not at all random, that
% exactly 50% of the elements in the first half are 1, that is, IF the
% sequence was truly random. Worse, try the same thing for ANY number of
% perfect riffle shuffles.
S = 20;
B = B0;
for ns = 1:S
B = B(riffleperm);
disp([ns,sum(B(1:N))])
end
1 250000
2 250000
3 250000
4 250000
5 250000
6 249984
7 250000
8 250000
9 250000
10 250000
11 250000
12 249996
13 250000
14 250016
15 250000
16 250002
17 250000
18 249998
19 249992
20 249996
% Surprisingly, a fair number of those shuffles end up with exactly 50% of
% the elements 1 in each half o the sequence. Perhaps mor interesting yet,
% is to look at the number of elements in either half after EXACTLY 180
% such perfect riffle shuffles.
S = 180;
B = B0;
for ns = 1:S
B = B(riffleperm);
end
sum(B(1:N))
ans =
0
%% Yes, after 180 perfect riffle shuffles, that sequence of 1e6 bits
% really has returned to its original order. Had you started with
% B = 1:1e6;
% you would indeed find the final sequence in B is perfectly 1:1e6 again.
% 180 perfect riffle shuffles are sufficient there to completely and
% exactly restore the sequence of 1:1e6.
B = 1:2*N;
S = 180;
for ns = 1:S
B = B(riffleperm);
end
isequal(B,1:2*N)
ans =
logical
1
Now you should recognize something interesting. That is, after 1 perfect riffle shuffle, we will have identically the same sequence as 181 perfect riffle shuffles, since 180 shuffles returned us to the starting position. In fact, every 180 such perfect riffle shuffles are sufficient to return the sequence to the exact start point.
The most interesting question in all of this might be to show that for N = 5e5, thus a total sequence of 1e6 such elements, that 180 is the magic number of perfect riffle shuffles to restore the original sequence, when N = 5e5.
That probably would get more deeply into group theory then we want to go here.
4 Comments
John D'Errico
on 17 Jul 2020
Edited: John D'Errico
on 17 Jul 2020
Your belief structure is not that relevant.
First, you need to accept that every perfect shuffle involves exactly the same resequencing of the cards. That is why the indexing solution I showed really is a valid shuffling sche,e. It is created by performing the shuffle on the numbers 1:2*N, Subsequent shuffles are the exact same permutation.
In fact, we can view a shuffle as a matrix multiply. This is why the corresponding matrix is called a permutation matrix. (A standard idea in mathematics.)
Consider a deck of size 8, so as you had it, N = 4.
N = 4;
% riffleperm has the property that B(riffleperm) will be a perfect riffle
% shuffle.
riffleperm = zeros(1,2*N);
riffleperm(1:2:end) = 1:N;
riffleperm(2:2:end) = N+1:2*N;
If you look at the sequence created, we see:
riffleperm
riffleperm =
1 5 2 6 3 7 4 8
In fact, that is exactly the out-shuffle you were using, but because I use a vector index, this is a faster way of doing it.
Now consider how it applies to the vector 1:8. According to the wikipedia page, on an 8 card deck, exactly 3 out-shuffles will restore the deck. We can find that same informatino from the OEIS - the online encyclopedia of integer sequences. Does it really work?
B = 1:8;
% after one shuffle:
B = B(riffleperm)
B =
1 5 2 6 3 7 4 8
% The second shuffle:
B = B(riffleperm)
B =
1 3 5 7 2 4 6 8
% and the third shuffle restores the sequence
B = B(riffleperm)
B =
1 2 3 4 5 6 7 8
So Wikipedia and the OEIS got it right. 3 perfect shuffles restore the sequence.
By the way, we can also view this as a matrix multiply. Again, they call it a permutation matrix for a good reason. In fact though, this would probably be a really inefficient way to do the permutation, as a matrix multiply, even if we used sparse matrices.
P = eye(8);
>> P = P(:,riffleperm)
P =
1 0 0 0 0 0 0 0
0 0 1 0 0 0 0 0
0 0 0 0 1 0 0 0
0 0 0 0 0 0 1 0
0 1 0 0 0 0 0 0
0 0 0 1 0 0 0 0
0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 1
>> B*P
ans =
1 5 2 6 3 7 4 8
>> B*P^2
ans =
1 3 5 7 2 4 6 8
>> B*P^3
ans =
1 2 3 4 5 6 7 8
Yes, there it is again. 3 applications of P to the original sequence restored the sequence. Mathematics is a pretty thing some days. A nice way to visualise this is to see that for out-shuffle permutation, we have P^3 equals the 8x8 identiy matrix.
P^3
ans =
1 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0
0 0 1 0 0 0 0 0
0 0 0 1 0 0 0 0
0 0 0 0 1 0 0 0
0 0 0 0 0 1 0 0
0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 1
Note that not ALL perfect shuffles of an 8 card deck restore the vector in exactly 3 shuffles. For example, I'll generate a random permutation.
riffleperm = randperm(8)
riffleperm =
3 7 1 5 8 4 2 6
B
B =
1 2 3 4 5 6 7 8
P = eye(8);
P = P(:,riffleperm);
B*P
ans =
3 7 1 5 8 4 2 6
B*P^2
ans =
1 2 3 8 6 5 7 4
B*P^3
ans =
3 7 1 6 4 8 2 5
B*P^4
ans =
1 2 3 4 5 6 7 8
So with this other specific shuffle, 4 perfect shuffles of that form are sufficient to restore the 8 card deck. A 5th shuffle is equivalent to what you get after exactly one shuffle. We can show this by understanding that P^4 power is the identity matrix. Therefore, we should see that P^5=P^4*P, so a 5th such shuffle of this form will be equivalent to the first shuffle.
P^4
ans =
1 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0
0 0 1 0 0 0 0 0
0 0 0 1 0 0 0 0
0 0 0 0 1 0 0 0
0 0 0 0 0 1 0 0
0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 1
I hope you understand how these things work, that is permutations and how to view a perfect shuffle as a permutation, and finally how to incorporate that into a permutation matrix. The nice thing about the permutation matrix form is it is easy to understand how a permutation can be viewed as a matrix multily, and how powers of the permutation matrix then result in a repeated shuffle. It all holds together rather neatly. As I said before, we could even get more deeply into the group theory concepts of such things.
The same ideas apply to larger "decks", even a deck with 1e6 total elements in it. As it turns out, with 1e6 total elements, that specific shuffle does indeed repeat after 180 permutations. Of course, I don't want to raise a 1e6 square matrix to the 180'th power. but the vector index form is sufficient, and I imagine it would be faster. :)
N = 5e5;
riffleperm = zeros(1,2*N);
riffleperm(1:2:end) = 1:N;
riffleperm(2:2:end) = N+1:2*N;
B0 = 1:2*N;
B = B0;
for i = 1:180
B = B(riffleperm);
end
isequal(B,B0)
ans =
logical
1
Again MATLAB sees the result after exactly 180 shuffles are sufficient to recover the original sequence. The 200th shuffle would then be the same sequence as the 20th shuffle.
B20 = B0;
for i = 1:20
B20 = B20(riffleperm);
end
B200 = B0;
for i = 1:200
B200 = B200(riffleperm);
end
isequal(B20,B200)
ans =
logical
1
Again, this is how it works. A belief structure is not necessary, because you can always prove it to be fact. That is the beauty of mathematics.
Now, how can we show that 180 is the recovery "time" for a shuffle of this form on 1e6 total elements? Again, the wiki page suggests this is based on the multiplicative order of 2, modulo 2*N-1. 2*N-1 is used here because your shuffle involved 2*N cards. Effectively, this asks us to compute the smallest value of X, such that
mod(2^X,1e6-1) == 1
Be careful here, as large powers of 2 will rapidly exceed the dynamic range of double precision arithmetic. 2^180 is far beyond the limit, which is 2^53-1 to represent an integer exactly. However MATLAB offers the powermod function, which allows us to compute the necessary remainder as:
powermod(2,180,1e6-1)
ans =
1
Indeed, we see the remainder is 1. 180 is the smallest such power of 2 that does what we need.
Going back to the 8 card deck, all we need to know is:
mod(2^3,8 - 1)
ans =
1
So for that size deck, exactly 3 out-shuffles really should recover the original sequence.
It takes only a wee bit of work to show, for example, that given a deck of size
D = 12345678;
powermod(2,10136,D-1)
ans =
1
Therefore, given that size deck, exactly 10136 perfect out-shuffles will resoire the original deck to perfect order.
Again, belief is not necessary, just mathematics. And with the correct mathematics, you can prove what you need.
More Answers (2)
Steven Lord
on 12 Jul 2020
function out = shuffle(in)
out = in; % Make out the correct size and type first
out(1:2:end) = in(1:end/2);
out(2:2:end) = in((end/2)+1:end);
end
This is for an out shuffle on a deck with an even number of cards. Modifying it to handle an in shuffle is an exercise for the reader, as it modifying this to work for a deck with an odd number of cards (if required by the terms of your assignment, like if you wanted to shuffle the thirteen cards of an individual suit rather than all fifty-two cards of a standard poker deck.)
x = 1:10;
for numShuffles = 1:6
x = shuffle(x);
fprintf('The deck after shuffle %d is %s.\n', numShuffles, mat2str(x))
end
0 Comments
See Also
Categories
Find more on Particle & Nuclear Physics in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!