MATLAB Answers

When does randomness enter a shuffled "deck" of coin tosses.

90 views (last 30 days)
george simpson
george simpson on 12 Jul 2020 at 13:32
Commented: george simpson on 17 Jul 2020 at 11:23
%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

David Hill
David Hill on 12 Jul 2020 at 15:13
I do not understand your question completely. Why not just:
A=[ones(1,5e5),-ones(1,5e5)];
A=A(randperm(1e6));%for each shuffle
John D'Errico
John D'Errico on 12 Jul 2020 at 17:45
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.
george simpson
george simpson on 12 Jul 2020 at 22:00
david hill and john D'errico ;
Thank you. I was hoping to
riffle shuffle my deck of 1e6 clean cointosess of +/-1's ,and then ultimately apply "runstest" for an exercise.
I have manually cranked to >200 shuffles, and each show a different Random Walk.All of my cointosses will have an equal # of heads and tails, but different RW's. #200 RW resembles a sine curve.
The help ii'm looking for is how to loop the operation so that I can get to a very large number of shuffles efficiently/quickly. At 200 shuffle i need 10 s of machine time.
George

Sign in to comment.

Accepted Answer

John D'Errico
John D'Errico on 12 Jul 2020 at 22:30
Edited: John D'Errico on 12 Jul 2020 at 22:39
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

Show 1 older comment
John D'Errico
John D'Errico on 13 Jul 2020 at 13:59
There are also subtly alternative ways to shuffle a deck. But all of them can be turned into a shuffle as I did, by resequencing an index vector, and then using that index vector to perform the shuffle.
At a glance, the riffle shuffle you are using looks like an out shuffle. In those pages, you will see the number of perfect shuffles of that sort to return a deck with an even number of cards to its original order. In fact, they are even tablulated in the OEIS. Thus, for the outshuffle, apprarently it is sequence A002326 in the OEIS.
And again, once you have reached the number of shuffles to return the deck to order and the deck is in its original sequence, then futher shuffles act as if you started all over again at the very beginning.
An important thing to recognize is that for the problem at hand, since only 180 shuffles are used to restore the deck to perfect order, that there are only 180 distinct possible permutations of this deck, when using out shuffles. Of course, there are a vast number of possible permutations. For example, the number of ways we could shuffle those 1e6 bits is literally immense. With N = 5e5, we would have
factorial(2*N)/factorial(N)^2
We cannot compute that in MATLAB using double precision. Even using symbolic computations, it would be a number with a VAST number of digits. We can use the gamma function to compute the necessary power of 10 as:
>> (gammaln(2*N + 1) - 2*gammaln(N+1))/log(10)
ans =
300966.691648236
So a number on the order of 10^300966. there are roughly 301 THOUSAND decimal digits in the number of ways we could permute those bits. That is one BIG number.
However, we can never realize more than 180 possible such permutations using the riffle shuffles you have indicated. And that is effectively one SMALL number in context. So I don't think you would ever really call any permutation you will get here a truly random sequence. :) The only question is, how non-random is it, or how random does it look.
george simpson
george simpson on 13 Jul 2020 at 19:27
John D'errico;
THanks for sharing your continued thoughts on the sbject of shuffling leading toward whatever may be meant by a perfect shuffle. There is some value in it, but it is not in the direction I would follow. I am an experimentalist, and tend to avoid theorizing as you are doing, So, if you want pick up the ball and run with it, I have no objections, nor could I effectively support you in all of your arguments.
I'm not convinced that any continued shuffle starting from a "clean" deck of 5e5 heads and 5e5 tails will cycle as you say: "clean" --> almost random--> to almost "clean". That is why I am playing this exercise. To see what I could see.
The almost clean may involve formation of longer, and longer Runs among all the possibilites occurring .
I have made better that 200 of these shuffles so far . When I start analyzing the results, I'll be on the look out for almost random. Maybe 180 will do the trick. I'll try to let you know what I find.
George
John D'Errico
John D'Errico on 17 Jul 2020 at 2:01
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.

Sign in to comment.

More Answers (2)

Steven Lord
Steven Lord on 12 Jul 2020 at 17:26
You're trying to simulate a perfect riffle shuffle?
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

Sign in to comment.


george simpson
george simpson on 12 Jul 2020 at 22:18
Steven lord;
Thank you. I will study your command, But ...
What i want to do is "shuffle" a deck of 1e6 +/-1: cointosses. My shuffle command is adequate,
to start with a clean deck, and study the changes in random walks from the shuffled coin tosses.
I would like to loop the calculations of shuffle operation so that I may do 1e3 ahuffles in succession.
I can do it by hand, but slowly and chews up machine time,
Maybe, I can find a number of shuffles that will randomize the random walk. If 7 shuffledswork for cards,
will something work for RW's
George

  6 Comments

Show 3 older comments
Steven Lord
Steven Lord on 16 Jul 2020 at 20:54
I suspect your question is along the line of “who is george simpson “, “what is he up to’, and how to get him to solve his own problem” ,’does he have a license to drive his Malab”,….
Do try accepting...
No. I want to know what problem you're trying to solve to understand if there's a more effective tool already in MATLAB that we can suggest to solve the problem you're trying to solve.
If you ask the question "How do I dig a hole?" the answer will differ depending on how big a hole you want. If you're planting flower bulbs, the answer is "Use a garden trowel." If you're putting in a tree, you probably want a shovel. If you're trying to lay foundation for a new house, the right tool for the job is probably an excavator.
But to answer the final question you asked, don't define new variables C2, C3, C4, etc. each time.
deck = [ones(1, 10), -ones(1, 10)];
numShuffles = 16;
history = repmat(deck, numShuffles+1, 1);
for whichDeckToProcess = 1:numShuffles
currentDeck = history(whichDeckToProcess, :);
A = currentDeck(1, 1:10);
B = currentDeck(1, 11:20);
C = reshape([A; B], 1, 20);
history(whichDeckToProcess+1, :) = C;
end
george simpson
george simpson on 16 Jul 2020 at 21:03
Steven Lord:
Thank you I'll try to make that wotk.
George
george simpson
george simpson on 17 Jul 2020 at 11:23
July 17,2020
John D’Errico;
Thanks for the extremely thoughtful reply about your solution to the return of a 1e6 deck, along with relevancy of belief systems.
I promise, I will study your arguments,eventually something will sink in. But mostly ,
I will try to learn from your Matlab commands, how you accomplish what you do.
I suspect that something got missed in the correspondence. I haven’t been asking for “HELP” in solving “ when a 1e6 deck is shuffed until it’s random”, nor, “how to make a perfect shuffle”,but how to make a “for loop”!
I “believe”, what you say is informative, and relates to question of what it takes to return to original shuffle. I believe it is publishable. But, It just doesn’t do it for me. I am not a mathematician.
By the way,I’m convinced ,now, that the shuffle of even 1e6 cards will eventually get back to shuffle #1, mostly, by my shuffle of a ten card deck.
.
What I am trying to do is show by “ experiment” that the shuffled clean deck will
become random.
But , I didn’t ask for help in solving that, that is what I’m currently working on. I asked about making a “for loop” using my shuffle to examine up to 1e6 shuffles.
{ Steven Lord gave me one, that I am testing, but if you could supply another, great. Please do so.}
Your arguments , while publishable, aren’t what I’m looking for.
By the way, In my shuffles, #180 isn’t a return , but is a “random” by runstest standards. There seems to be an early return around 5 or 6, that I’ll study further, but for me the return is moot...
George

Sign in to comment.

Products