Within any $k$ cards there are $\binom{k}{2} = k(k-1)/2$ pairs. Thus, in a "shoe" that has $m$ distinct cards with $k_1$ of one kind, $k_2$ of another, and so on, there are $$\binom{k_1}{2} + \binom{k_2}{2} + \cdots + \binom{k_m}{2}$$ pairs of like cards and there are $$\binom{k_1+k_2+\cdots+k_m}{2}$$ pairs of cards, like or unlike, altogether. When drawing two cards at random, each of these pairs has an equal chance. Therefore:
The chance of drawing a pair of like cards is $$\frac{\binom{k_1}{2} + \binom{k_2}{2} + \cdots + \binom{k_m}{2}}{\binom{k_1+k_2+\cdots+k_m}{2}}.$$
After drawing two cards--whether or not they are like--update the $k_i$ accordingly and repeat.
For example, in four decks of playing cards there are (presumably) $m=13$ distinct cards and $16=k_1=k_2=\cdots=k_m$ of each. At the outset, then, the chance of drawing a pair of like cards is
$$\frac{13\times\binom{16}{2}}{\binom{13\times 16}{2}} = \frac{1560}{21528} \approx 7.2464\%.$$
If the first two cards were a like pair, then one of the $k_i$ is reduced from $16$ to $14$ and the new chances on the next draw are
$$\frac{12\times\binom{16}{2} + \binom{14}{2}}{\binom{12\times 16\ +\ 14}{2}} = \frac{1531}{21115} \approx 7.2508\%.$$
Otherwise, if the first two cards were not a like pair, then two different $k_i$ are reduced from $16$ to $15$ and the new chances on the next draw are
$$\frac{11\times\binom{16}{2} + 2\times \binom{15}{2}}{\binom{11\times 16\ +\ 2\times 15}{2}} = \frac{1530}{21115} \approx 7.2460\%.$$
If you need to keep going with this calculation as cards are drawn, evidently you must keep track of all the $k_i$--there is no shortcut.
The algorithm is particularly simple, though, because when you draw a card from a group of $k_i$ like cards, $\binom{k_i}{2}$ decreases by $k_i-1.$
The following algorithm maintains a data structure containing (a) a vector of the $k_i,$ (b) the numerator and denominator of the fraction, and (c) the sum of the $k_i$ (let's call this $n$).
At the beginning you have to use the binomial formula above to compute the numerator $a$ and denominator $b;$ the fraction $a/b$ is the chance of drawing a pair.
Before drawing two cards the chance of a like pair is $a/b.$ Each time a card is drawn from group $i,$ decrease $a$ by $k_i-1,$ decrease $b$ by $n-1,$ decrease $n$ by $1,$ and finally decrease $k_i$ by $1.$
An R
implementation is shown below to illustrate. Here is part of its output:
a = 1560 ; b = 21528 ; Chance of a pair is 0.07246 Drew 7 and 6
a = 1530 ; b = 21115 ; Chance of a pair is 0.07246 Drew 6 and 4
a = 1501 ; b = 20706 ; Chance of a pair is 0.07249 Drew 6 and 6 (pair)
a = 1476 ; b = 20301 ; Chance of a pair is 0.07271 Drew 3 and Q
a = 1446 ; b = 19900 ; Chance of a pair is 0.07266 Drew K and K (pair)
... (98 missing lines)
a = 0 ; b = 1 ; Chance of a pair is 0.00000 Drew 3 and 2
Here is code that carries out the algorithm and an example of its use.
deck <- rep(c("A", 2:10, c("J", "Q", "K")), 4)
shoe <- rep(deck, 4)
#
# Initialize the data structure `X`.
#
k <- table(shoe)
n <- sum(k)
a <- sum(k * (k-1) / 2)
b <- n * (n-1) / 2
X <- list(k=k, n=n, a=a, b=b)
#
# Here is how to update the data structure when a card `i` is drawn.
#
update <- function(i, X) {
i <- as.character(i)
X$a <- X$a - (X$k[i]-1)
X$b <- X$b - (X$n-1)
X$n <- X$n - 1
X$k[i] <- X$k[i] - 1
return(X)
}
#
# Illustrate.
#
set.seed(17)
draws <- sample(shoe, 208) # Shuffle the cards
for (j in 1:(length(draws)/2)) { # Draw two at a time
with(X, cat("a =", a, "; b =", b, "; Chance of a pair is", sprintf("%.5f", a/b)))
X <- update(draws[2*j-1], X)
X <- update(draws[2*j], X)
with(X, cat("\tDrew", draws[2*j-1], "and", draws[2*j],
ifelse(draws[2*j-1]==draws[2*j], "(pair)", ""), "\n"))
}