I have multiple stochastic processes that can have 2 different realisations. I think of them as coins that can take the values heads or tails. For the example I assume that there are 3 coins. I know the 3 probabilities of showing heads that these coins have, but I don't know which belongs to which coin. For example I know that the probabilities are 60%, 50%, and 40% and each occurs exactly once i.e. if Coin A has a probability of 60% for heads, Coin B and C cannot have 60%.

I observe realisations for each coin over multiple rounds. With 3 coins there are 3! = 6 different ways the probabilities could be assigned to the coins. For each of these 6 possibilities I want to know probability that it is the true one given the data on the realisations of the 3 coins that I have. For example I want a result like with 20% probability Coin A has 60% probability of heads, Coin B has 50% and Coin C has 40% and with 15% Probability Coin A has 50% probability of heads, Coin B has 60% and Coin C has 40% and so on.

How can I calculate these probabilities?

  • $\begingroup$ interesting question. It's not hard to come up with the most likely ordering, but deriving probabilities seems harder. Actually, with a computer, it's not so hard I think $\endgroup$ – Cam.Davidson.Pilon Aug 31 '18 at 19:16
  • $\begingroup$ The nature of your situation is unclear: I cannot see the sense in which "2 different realizations," "3 coins," and "multiple rounds" all describe the same experiment. Could you perhaps display and explain a small example? $\endgroup$ – whuber Sep 3 '18 at 13:46
  • $\begingroup$ To clarify: the two different possible realisations are heads or tails if we thing of the processes as coins. I have 3 biased coins that I flip at the same time and observe the results. After I have observed the results I flip all 3 coins again and observe the results. This is what I mean by multiple rounds. So if I have 4 rounds and 3 coins for example that would mean that I flip each of the 3 coins 4 times and observe heads or tails $\endgroup$ – umbal Sep 5 '18 at 8:12

The idea is to iterate over all possible permutations, compute the posterior probability (likelihood * prior), and then normalize.

Here's a Python program to compute the posterior probabilities:

from math import factorial
from itertools import permutations

# these are our fixed probabilities
P = [0.4, 0.5, 0.6]

# need a prior
PRIOR = 1.0/3.0

# here's the observed data: 5 "rounds", and the number of successes of each coin.
observed_data = (5, (3,2,1))

# some helper functions
def C(n, r):
    return factorial(n) / factorial(r) / factorial(n-r)

def likelihood(observed_data, permutation):
    permutation is a permutation of (0,1,2)


    N, successes = observed_data
    prob = 1.0
    for observed_i, proposed_i in enumerate(permutation):
        p = P[proposed_i]
        r = successes[observed_i]
        prob *= C(N, r) * p**r * (1 - p)**(N - r)
    return prob 

def unnormalized_posterior_probability(observed_data, permutation):
    observed_data is a tuple, first element is the number of rounds
    the second element is a list representing the count of successes
    ex: (5, [2, 0, 3])

    return likelihood(observed_data, permutation) * PRIOR 

# compute the posterior probabilities
running_sum = 0
posterior_probabilities = {}

for possible_order in permutations((0,1,2)):
    posterior_probabilities[possible_order] = unnormalized_posterior_probability(observed_data, possible_order)
    running_sum += posterior_probabilities[possible_order]

for possible_order in permutations((0,1,2)):
    posterior_probabilities[possible_order] = posterior_probabilities[possible_order] / running_sum

The output looks like:

{(0, 1, 2): 0.063,
 (0, 2, 1): 0.094,
 (1, 0, 2): 0.094,
 (1, 2, 0): 0.213,
 (2, 0, 1): 0.213,
 (2, 1, 0): 0.320}

Which suggests that (2, 1, 0) is the most probable order (which makes sense, given P and the observed data.

  • $\begingroup$ Thank you for the code, this is an improvement over my own approach.I think it is still not the ideal solution because it does not use all information. The issue I see is that in the "def likelihood:" part you get the likelihood of a permutation by calculating individual likelihoods for each of the 3 coins and then multiply them. I think that relies on a independence assumption that we can't make. When we look at the data and ask how likely is it that a coin with 0.4 probability generated the data we need to account for the fact otherwise it must have been one with 0.5 or 0.6 probability. $\endgroup$ – umbal Sep 5 '18 at 10:16
  • $\begingroup$ > we need to account for the fact otherwise it must have been one with 0.5 or 0.6 probability. hm, that part is done in the normalization in the for loops, no? $\endgroup$ – Cam.Davidson.Pilon Sep 5 '18 at 10:54
  • $\begingroup$ My intuition is that taking the fact that the only possible probabilities for the coins are 0.4,0.5 and 0.6 into account in the step where you calculate the probability of a permutation is not necessarily the same as calculating the probabilities for each coin individually, multiply them up and then normalise later. I might be wrong about this but i think with Bayesian updating their should be a way to do this without normalising. $\endgroup$ – umbal Sep 5 '18 at 16:51
  • $\begingroup$ My interpretation: the information is being used, but it's implicit. Going back to Bayes: $P(A=a|B=b) = P(B=b | A=a) P(A=a) / P(B=b)$ - note that our likelihood term, $P(B=b | A=a)$ we are using that $A=a$, and implicitly that $A \ne a', A \ne a''$, etc. So we are using the information in this sense.Similarly, in our example, in $P(B=b)$ we are using that there are only three possibilities for $A$, again implicitly: $P(B=b) = P(B=b | A=a) P(A=a) + P(B=b | A=a') P(A=a') + P(B=b | A=a'') P(A=a'')$. $\endgroup$ – Cam.Davidson.Pilon Sep 5 '18 at 19:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.