## Dependent urns problem

This problem has been suggested by my professor.
Maybe it is trivial but I am not able to solve it. Thank you in advance for your help!

At time t=0, there are 4 urns: A, B, C, D. Every urn contains a positive number of red and black balls.
Urn A contains b_A0>0 black balls and r_A0>0 red balls. Urn B b_B0>0 and r_B0>0, etc.
Let n_A0=b_A0+r_A0, n_B0=b_B0+r_B0, etc...

The four urn are connected as follows

A----B----C----D

We focus our attention on urns B and C
Urn B has two neighbors: A and C
Urn C has two neighbors: B and D

At time t=1, the two urns are sampled according to the following rule:
Urn B
- a ball is sampled at random from the urn and the probability that this ball is black is not b_B0/n_B0, but (b_B0+b_A0+b_C0)/(n_A0+n_B0+n_C0). In other words, the probability of picking a black ball from urn B depends on the total number of black balls in its neighborhood. It is like A, B and C form one big urn, containing all their balls.
An analogous reasoning is valid for red balls.
- once we have discovered the color of the sampled ball, we replace it into urn B and we add one extra ball of the same color (as in a Polya urn) to B.
- if X_t is the color of the sampled ball at time t and we assume that X_t=1 means "black ball sampled", we have that X_t is a Bernoulli variable with parameter (b_Bt-1+b_At-1+b_Ct-1)/(n_At-1+n_Bt-1+n_Ct-1).
We have that b_Bt=b_Bt-1+1if X_t=1 and b_Bt=b_Bt-1 if X_t=0.

Urn C
Always at time t=1, we also sample urn C (according to the same mechanism with reinforcement). In this case we consider the variable Y_t, which indicates the sampled ball from urn C at time t.
Even Y_t is a Bernoulli variable with parameter (b_Bt-1+b_Ct-1+b_Dt-1)/(n_Bt-1+n_Ct-1+n_Dt-1).

My problem is the following: it is obvious that the two urns are not independent. In fact the evolution of black (or red) balls in urn A also depends on the composition of urn B and viceversa.
What is the joint probability of X_t and Y_t, i.e. P(X_t=x;Y_t=y) where x,y=0,1?
Notice that for t=0, P(X_0=x;Y_0=y)=P(X_0=x)P(Y_0=y), i.e. the product of two Bernoulli variables.
But this is not true for t>=1.

Thanks!!!