Anomalous probability result for computed single pair probability when dealing seven card poker hands

Question

BACKGROUND: Today I thought I'd begin a small project building a poker simulator. The first task I set was to deal cards from a shuffled deck, and check the various numerically generated probabilities against accepted values. The first such probability I checked was the single pair probability--that is, generating (numerically) the probability of being dealt a single pair, given as inputs the number of cards dealt and the number of hands dealt, where each hand is dealt from a separate shuffled deck. Cards are dealt from the top of the deck. Below I show the beginning of that program. I first tested numerically generated single pair probability for five card hands. The computed value comes to within a tenth of a percent of the accepted single pair probability for five card hands (but always high by about a tenth of a percent): https://en.wikipedia.org/wiki/Poker_probability

However, when I test the numerically generated single pair probability for seven card hands, I find that I am off by 4% to 5% from the accepted value (e.g., typical computed value = 0.47828; accepted value as per above = 0.438). I've run the numerical experiments up to ten million hands dealt. The computed single pair probability for seven card hands is stable, and remains off by 4% to 5% from the accepted value. It's not clear why this is the case.

QUESTION: Why is this the case? I suspect that my code is not taking something into account, but I cannot detect what. Python code follows . . .

NOTE: Issue 31381901 is similar to this one. But in the below code the issue of double counting is accounted for by converting the dealt hand to a set, which will eliminate duplicate values, thus reducing the size of the set (in the case of 7 card hands) from 7 to 6. That reduction indicates a single pair. If three-of-a-kind is present, then the size of the set would be 5, since two of the three cards in the three-of-a-kind would be eliminated by the set conversion.

from random import shuffle

def make_deck():
    '''Make a 52 card deck of cards.  First symbol
    is the value, second symbol is the suit.  Concatenate
    both symbols together.
    Input:  None
    Output:  List
    '''
    value = ['A', '2', '3', '4', '5', '6', '7', '8', '9', 'T', 'J', 'Q', 'K']
    suit = ['C','H','S','D']
    deck = [j+i for j in value for i in suit]
    return deck


def shuffle_deck(deck, times_to_shuffle=7):
    '''Shuffle a deck of cards produced by make_deck().
    Default:  7 times.
    Input:  list, int
    Output:  list (modified in-place)
    '''
    for n in range(times_to_shuffle):
        shuffle(deck)


def test_for_single_pair(hand, cards_per_hand):
    '''Tests for presence of a single pair in 
    a dealt hand by converting the hand to a set.
    The set representation of a hand with a single
    pair will have one less member than the original
    hand.
    Input: list, int
    Output:  int
    '''
    hand_values_lst = [card[0] for card in hand]
    hand_values_set = set(hand_values_lst)
    set_size = len(hand_values_set)

    if set_size == (cards_per_hand - 1):
        return 1
    else:
        return 0


def deal_series_of_hands(num_hands,cards_per_hand):
    '''Deals a series of hands of cards and tests
    for single pairs in each hand.  Creates a deck
    of 52 cards, then begins dealing loop.  Shuffles
    deck thoroughly after each hand is dealt.
    Captures a list of the dealt hands that conform
    to the spec (i.e., that contain one pair each),
    for later debugging purposes
    Input:  int, int
    Output:  int, int, list
    '''

    deck = make_deck()
    single_pair_count = 0
    hand_capture = []

    for m in range(num_hands):
        shuffle_deck(deck)
        hand = deck[0:cards_per_hand]    #first cards dealt from the deck
        pair_count = test_for_single_pair(hand, cards_per_hand)
        if pair_count == 1:
            single_pair_count += pair_count
            hand_capture.append(hand)

    return (single_pair_count, num_hands, hand_capture)



cards_per_hand = 7   #User input parameter
num_hands = 50000    #user input parameter

single_pair_count, num_hands_dealt, hand_capture = deal_series_of_hands(num_hands, cards_per_hand)
single_pair_probability = single_pair_count/ num_hands_dealt
single_pair_str = 'Single pair probability (%d card deal; poker hands): '%(cards_per_hand)
print(single_pair_str, single_pair_probability)

Answer 1

If the hand contains a single pair but also contains a higher-value unit such as a straight or a flush, your code still counts that as a pair, where the probability article does not.