I spent most of this week working on joint probability spaces, and have uploaded an incomplete pull request, #14764. This PR implements a basic structure of the joint probability classes. Like I said, the PR is still incomplete and requires a lot of work before it can be merged, but hopefully, I should be able to complete the implementation on time. Right now, the work done primarily focusses on creating Joint probability spaces and the underlying domain and density out of independent random variables. Example:

>>> from sympy.stats.joint_rv import *
>>> from sympy.stats import Geometric, Poisson
>>> from sympy import S
>>> from sympy.abc import x, y
>>> X, Y = Geometric('X', S(1)/2), Poisson('Y', 4)
>>> Z = Joint('Z', (X, Y))
>>> density(Z)(x, y)
2**(-x + 1)*4**y*exp(-4)/(2*factorial(y))

For the next week, I plan to implement joint spaces so as to include user defined distributions, for which individual random variables have not been declared earlier.