Let's have an urn with $m_i= 2$ balls for $c=3$ different colours with different weights $w_i$.
$n=2$ balls are taken randomly but the probability of sampling a particular coloured ball is proportional to its weight. This is a biased urn problem and the probability of taking a given sample follows the Multivariate Fisher's Noncentral Hypergeometric Distribution with the following probability mass function ($S$ is the support of the PMF):
$$ \frac{1}{P_0}\prod_{i=1}^{c} \binom{m_i}{x_i}\omega_i^{x_i}\\ P_0 = \sum_{(y_0,\ldots,y_c)\in \mathrm{S}}\prod_{i=1}^{c} \binom{m_i}{y_i}\omega_i^{y_i}\\ \mathrm{S} = \left\{ \mathbf{x} \in \mathbb{Z}_{0+}^c \, : \, \sum_{i=1}^{c} x_i = n \right\} $$
I would like to calculate the probability of the realization $x = \{1, 0, 1\}$.
The R
's package BiasedUrn
exposes the Multivariate Fisher's Noncentral Hypergeometric Distribution. Here a reproducible example:
library(BiasedUrn)
set.seed(123)
n = 2 # balls taken
c = 3 # different colors
m = rep(2, c) # nuber of balls in the urn for each colour
w = runif(n = c)
w = w/sum(w) #weights of the balls (biasd urn)
# a realization
x = c(1, 0, 1)
prob = dMWNCHypergeo(x, m, n, odds = w)
The calculated probability is: $prob = 0.1209264$.
when I try to calculate it with the formula above I get different results:
numerator = 1.0
for (i in 1:c) {
numerator = numerator * choose(m[i], x[i]) * (w[i] ^ x[i])
}
S = combn(x=rep(c(1:c),n), m = n, tabulate, nbins = c)
S = S[, !duplicated(t(S))]
denominator = 0.
for (y in S) {
prod = 1
for (i in 1:c) {
prod = prod * choose(m[i], y) * (w[i]^y)
}
denominator = denominator + prod
}
prob_2 = numerator / denominator
and $prob\_2 = 0.1329621$. So, it's wrong!
While I am rather confident that the numerator is calculated in the right way, my instinct tells me that's something wrong in the denominator and probably I cannot identify in the right way the support of the PMF. I defined $S$ as the 6 possible (distinct) outcomes: $\{1,1,0\}$, $\{1,0,1\}$, $\{2,0,0\}$, $\{0,1,1\}$, $\{0,2,0\}$, $\{0,0,2\}$.
So, please, could someone help me find what I am doing wrong here?
I would like to generalize the problem and implement the Multivariate Fisher's Noncentral Hypergeometric Distribution in stan
to make inference on $\omega_i$ (an old question of mine) and so I need help to understand what's I am doing wrong here.