How to modify the Birthday Problem to answer the probability of 3, 4, 5 or more people in a room having the same birthday

Question

The Birthday Problem asks how many people need to be in a room for there to be greater than a 50% chance that at least two people have the same birthday. I've used the Python code below to solve for the probability of a least two people having the same birthday. I would like to modify the code so that I can find the probability of 3 or more, 4 or more (and so on) people having the same birthday.

I tried changing the variable c to 3 below but it didn't work

import random
import numpy as np
from math import factorial

round = 3    # precision of random numbers
n=10**round  # Number of possible permutations of the random numbers
k=50       # Number of samples in the analysis
c=2          # a constant represents every possible pairs (2) evaluated

# This is an empirical calculation of duplicated values
rvs = np.round(np.random.uniform(0,1,k), round)
unique=np.size(np.unique(rvs))
dups = (k-unique)


possible_combinations=n*(n-1)/2
print('n=', n)
print('k=', k)
print('Possible Pair Combinations = ', "{:,.0f}".format(possible_combinations))
print('Number of duplicate random values=', dups,'of ' +str(k) +' samples')

# Calculate probabilities (Method 1)
P_coeff=(k*(k-1)/c)
P_A = ((n-1)/n)**P_coeff
P_B =1-P_A

# Calculate probabilities (Method 2)
P_B2 = 1 - (factorial(n) / (factorial(n-k) * n**k))

print('P(A)=', "{:,.4f}".format(P_A))
print('P(B)=', "{:,.4f}".format(P_B)+ ' Method 1 - Probability of duplicated values')
print('P(B)=', "{:,.4f}".format(P_B2)+ ' Method 2 - Probability of duplicated values')

Answer 1

I finally found an answer to my questions at that uses a Poisson distribution: https://math.stackexchange.com/questions/25876/probability-of-3-people-in-a-room-of-30-having-the-same-birthday#:~:text=Then%20this%20approximation%20gives%20

I've written and tested Python code than answers the Birthday problem for 2, 3, 4 (etc) people having the same birthday:

# Solve for k people having the same birthday in a room with n people in it

from scipy.stats import poisson

# INPUTS
n=23      # Number of people in a room
k=2       # Number of people to evaluate for same birthday

# CALCULATIONS
k_=k-1     # k actuall need to be (k-1) assuming 1 person as the "seed"
b=365      # Days in a year
n_b=n/b    # mu for Poisson distribution

F_k = poisson.cdf(k_, n_b, loc=0)
F_k = F_k**b
P_k = 1 - F_k

print('Mu for Poisson CDF=',"{:,.4f}".format(n/b))
print('Poisson Prob of not having k birthdays=',"{:,.4f}".format(F_k))
print('Probability of k people in with birdays',"{:,.4f}".format(P_k))

Answers: Mu for Poisson CDF= 0.0630 Poisson Prob of not having k birthdays= 0.4988 Probability of k people in with birthdays 0.5012