We are going to show how we can estimate card probabilities by applying Monte Carlo Simulation and how we can solve them numerically in Python. The first thing that we need to do is to create a deck of 52 cards. Let’s start.

**How to Generate a Deck of Cards**

import itertools, random # make a deck of cards deck = list(itertools.product(['A', '2', '3', '4', '5', '6', '7', '8', '9', '10', 'J', 'Q', 'K'], ['Spade','Heart','Diamond','Club'])) deck

```
[('A', 'Spade'),
('A', 'Heart'),
('A', 'Diamond'),
('A', 'Club'),
('2', 'Spade'),
('2', 'Heart'),
('2', 'Diamond'),
('2', 'Club'),
('3', 'Spade'),
('3', 'Heart'),
('3', 'Diamond'),
('3', 'Club'),
('4', 'Spade'),
('4', 'Heart'),
('4', 'Diamond'),
('4', 'Club'),
('5', 'Spade'),
('5', 'Heart'),
...
]
```

**How to Generate Multiple Decks of Cards**

In some games, like Black Jack, it is possible to play with even 7 decks. Let’s see how we can generate 7 decks. We expect to see \(52\times 7 = 364\).

7deck = list(itertools.product([1,2,3,4,5,6,7], ['A', '2', '3', '4', '5', '6', '7', '8', '9', '10', 'J', 'Q', 'K'], ['Spade','Heart','Diamond','Club'])) len(d7eck)

`364`

**How to Shuffle the Deck**

# shuffle the cards random.shuffle(deck) deck[0:10]

```
[('6', 'Club'),
('8', 'Spade'),
('J', 'Heart'),
('10', 'Heart'),
('Q', 'Spade'),
('7', 'Diamond'),
('K', 'Diamond'),
('J', 'Club'),
('J', 'Diamond'),
('A', 'Club')]
```

**How to Sort the Deck **

For some probabilities where the order does not matter, a good trick is to sort the cards. The following commands can be helpful.

# sort the deck sorted(deck) # sort by face sorted(deck, key = lambda x: x[0]) # sort by suit sorted(deck, key = lambda x: x[1])

**How to Remove Cards from the Deck**

Depending on the Games and the problem that we need to solve, sometimes there is a need to remove from the Deck the cards which have already been served. The commands that we can use are the following:

# assume that card is a tuple like ('J', 'Diamond') deck.remove(card) # or we can use the pop command where we remove by the index deck.pop(0)

## Part 1: Estimate Card Probabilities with Monte Carlo Simulation

**Question 1: What is the probability that when two cards are drawn from a deck of cards without a replacement that both of them will be Ace? **

We can solve this probability explicitly:

\(\frac{4}{52}\times\frac{3}{51}=\frac{1}{221}=0.0045\)

Let’s say that we were not familiar with formulas or that the problem was more complicated. We could find an approximate probability with a Monte Carlo Simulation (10M Simulations).

N = 10000000 double_aces = 0 for hands in range(N): # shuffle the cards random.shuffle(deck) aces = [d[0] for d in deck[0:2]].count('A') if aces == 2: double_aces+=1 prob = double_aces/N prob

`0.0045214`

As we can see we found a good approximation. Good job. Let’s have a look at other problems.

**Question 2: What is the probability of two Aces in 5 card hand without replacement**.

Again, we can solve this explicitly by applying the formula of the hypergeometric distribution. We can use also Python for this calculation.

\(P = \frac{{4 \choose 2} { 48 \choose 3 }}{52\choose 5}\)

from scipy.stats import hypergeom # k is the number of required aces # M is the total number of the cards in a deck # n how many aces are in the deck # N how many aces cards we will get hypergeom.pmf(k=2,M=52, n=4, N=5)

`0.039929818081078594`

Let’s try to solve this numerically.

N = 10000000 double_aces = 0 for hands in range(N): # shuffle the cards random.shuffle(deck) aces = [d[0] for d in deck[0:5]].count('A') if aces == 2: # print(deck[0:2]) double_aces+=1 prob = double_aces/N prob

`0.0398805`

As we can see, again the estimated probability (0.0398805) was very close to the actual (0.0399298).

**Question 3: What is the probability of being dealt a flush (5 cards of all the same suit) from the first 5 cards in a deck? **

- The first card it does not matter what the suit is. Any of the suits can be drawn initially, as long as the next four cards are of the same suit as the original card.
- There are 13 of each suit in the deck, so after the first card is drawn, there are only 12 of that suit, then 11 left for the third card, 10 left for the fourth card, and 9 left for the final card. Also, there will be one less card total in the deck each time.

\(P(flush) = P( \textbf {2nd is same suit} ) P( \textbf {3rd is same suit} ) P( \textbf {4th is same suit} ) P( \textbf {5th is same suit} )\)

\(P(flush) =\frac{12}{51} \times \frac{11}{50} \times \frac{10}{49} \times \frac{9}{48} = \frac{33}{16660} = 0.00198 \)

Let’s solve this problem numerically. Here the trick will be search for 1 unique suit in the 5 cards.

N = 10000000 flushes = 0 for hands in range(N): # shuffle the cards random.shuffle(deck) flush = [d[1] for d in deck[0:5]] if len(set(flush))== 1: flushes+=1 prob = flushes/N prob

`0.0019823`

Again, we got a really good estimate!

**Question 4: ** **What is the probability of being dealt a royal flush from the first 5 cards in a deck?**

The number of al possible 5-card hand is \({52\choose 5} = 2,598,960\)

A royal flush is **ace, king, queen, jack, and ten of the same suit**. If we order the 5-card hand from highest card to lowest, the first card will be an ace. There are four possible suits for the ace. After that, the other four cards are completely determined. Thus, there are 4 possible royal flushes:

\(P(flush) = \frac{4}{ 2598960 } = 0.00000154 \)

Let’s see if the Monte Carlo works with such small numbers.

# royal flush N = 10000000 royal_flushes = 0 for hands in range(N): # shuffle the cards random.shuffle(deck) flush = [d[1] for d in deck[0:5]] face = [d[0] for d in deck[0:5]] if len(set(flush))== 1 and sorted(['A','J', 'Q', 'K', '10'])==sorted(face): royal_flushes+=1 prob = royal_flushes/N prob

`1.5e-06`

Awesome! Again we got a very good estimate!

## Part 2: Calculate Exact Card Probabilities Numerically

Above, we showed how we can calculate the Card Probabilities explicitly by applying mathematical formulas and how we can estimate them by applying Monte Carlo Simulation. Now, we will show how we can get the exact probability using Python. This is not always applicable but let’s try to solve the questions of Part 1.

The logic here is to generate all the possible combinations and then to calculate the ratio. Let’s see how we can get all the possible 5-hand cards of a 52-card deck.

# importing modules import itertools from itertools import combinations # make a deck of cards deck = list(itertools.product(['A', '2', '3', '4', '5', '6', '7', '8', '9', '10', 'J', 'Q', 'K'], ['Spade','Heart','Diamond','Club'])) # get the (52 5) combinations all_possible_by_5_combinations = list(combinations(deck,5)) len(all_possible_by_5_combinations)

`2598960`

As expected, \({52\choose 5} = 2,598,960\)

**Question 1: What is the probability that when two cards are drawn from a deck of cards without a replacement that both of them will be Ace? **

# importing modules import itertools from itertools import combinations # make a deck of cards deck = list(itertools.product(['A', '2', '3', '4', '5', '6', '7', '8', '9', '10', 'J', 'Q', 'K'], ['Spade','Heart','Diamond','Club'])) # get the (52 2) combinations all_possible_by_2_combinations = list(combinations(deck,2)) Aces = 0 for card in all_possible_by_2_combinations: if [d[0] for d in card].count('A') == 2: Aces+=1 prob = Aces / len(all_possible_by_2_combinations) prob

`0.004524886877828055`

**Question 2: What is the probability of two Aces in 5 card hand without replacement**.

# get the (52 5) combinations all_possible_by_5_combinations = list(combinations(deck,5)) Aces = 0 for card in all_possible_by_5_combinations: if [d[0] for d in card].count('A') == 2: Aces+=1 prob = Aces / len(all_possible_by_5_combinations) prob

`0.03992981808107859`

**Question 3: What is the probability of being dealt a flush (5 cards of all the same suit) from the first 5 cards in a deck? **

# get the (52 5) combinations all_possible_by_5_combinations = list(combinations(deck,5)) flushes = 0 for card in all_possible_by_5_combinations: flush = [d[1] for d in card] if len(set(flush))== 1: flushes+=1 prob = flushes / len(all_possible_by_5_combinations) prob

`0.0019807923169267707`

**Question 4: ** **What is the probability of being dealt a royal flush from the first 5 cards in a deck?**

# get the (52 5) combinations all_possible_by_5_combinations = list(combinations(deck,5)) royal_flushes = 0 for card in all_possible_by_5_combinations: flush = [d[1] for d in card] face = [d[0] for d in card] if len(set(flush))== 1 and sorted(['A','J', 'Q', 'K', '10'])==sorted(face): royal_flushes +=1 prob = royal_flushes / len(all_possible_by_5_combinations) prob

`1.5390771693292702e-06`

## Conclusion

Computing Power enables us to estimate and even calculate complicated probabilities without being necessary to be experts in Statistics and Probabilities. In this post, we provided some relatively simple examples. The same logic can be extended to more advanced and complicated problems in many fields, like simulating card games, lotteries, gambling games etc.

## 2 thoughts on “Estimate Probabilities of Card Games”