I've been thinking this through and so far I got this:
1) let's discuss a single variable transform. We have a random variable $X$ that transforms to $Y$ via formula $y = Y (x)$. The distribution for $X$ is $X \sim f_X (x)$, $Y \sim f_Y (y)$. The distributions are related as follows $$ f_X (x) \, \mathrm{d} x = f_Y (y) \, \mathrm{d} y $$ from which $f_Y (y)$ readily follows $$ f_Y (y) = \sum_i \left| \frac{\mathrm{d} x_i}{\mathrm{d} y} \right| f_X (Y_i^{-1}(y)) $$ where the complications are: map $x \to y$ might not be bijective, so there might be several branches of solutions, denoted by $i$, moreover, absolute value; the probabilities are added, regardless of the slope.
Now I actually came up with a somewhat more clever formula, that accounts for both of these complications and wraps up $f_Y (y)$ in a compact way:
$$ f_Y (y) = \int \limits_{-\infty}^\infty \mathrm{d} x \, f_X (x) \, \delta (y - Y (x)) $$
2) let's now approach a transformation, where several random variables are combined into another single random variable, i.e. $x, y \to z$ and $z = Z(x,y)$.
As of now, I am not aware of any formula of the form $f_Z (z) \, \mathrm{d} z = \cdots$, but my $\delta$-formula can be extended easily as follows $$ f_Z (z) = \int \mathrm{d} x \, \mathrm{d} y \, f_{X,Y} (x, y) \, \delta (z - Z (x, y)) $$
We can grasp the meaning in a two-step process: first find the accumulative function $F_Z (z)$ which just counts the total probability that $z < Z$ (we would have $\theta$ function instead of $\delta$ in the integral) and then we take the derivative of $F_Z (z)$ w.r.t. $z$ to obtain the formula above.
3) even more complicated case when several variables transform to several random variables $x, y \to z, w$ given by $z = Z (x,y)$, $w = W (x,y)$. My naive guess of the distribution for $z, w$ given distribution for $x, y$ is
$$ f_{Z, W} (z, w) = \int \mathrm{d} x \, \mathrm{d} y \, f_{X,Y} (x, y) \, \delta (z - Z (x, y)) \, \delta (w - W (x, y)) $$
My question is: is this correct? Can this be generalized even further, i.e. $$ f_{Y_1, \cdots Y_n} (y_1, \cdots, y_n) = \int \mathrm{d} x_1 \, \cdots \, \mathrm{d} x_m \, f_{X_1, \cdots, X_m} (x_1, \cdots, x_m) \times \\ \times \: \delta (y_1 - Y_1 (x_1, \cdots, x_m)) \, \cdots \, \delta (y_n - Y_n (x_1, \cdots, x_m)), $$ where $m \geq n$ ($m$ random variables transform into $n$ random variables)?
Thank you.
Yes, this is correct. However, whereas in the case $x,y\to z$ the delta distribution approach is useful in practice, in the case $x,y\to z,w$ it’s just a reformulation of the usual differential relationship that doesn’t add any practical value. We have
$$ f_{Z,W}(z,w)\mathrm dz\mathrm dw=f_{X,Y}(x,y)\mathrm dx\mathrm dy $$
and thus
$$ f_{Z,W}(z,w)=f_{X,Y}(x,y)\frac{\partial(x,y)}{\partial(z,w)}\;, $$
where $\frac{\partial(x,y)}{\partial(z,w)}$ is the Jacobian of the transformation. The delta distribution formulation follows directly by integration:
\begin{eqnarray} f_{Z,W}(z',w') &=& \iint\mathrm dz\mathrm dwf_{Z,W}(z,w)\delta(z'-z)\delta(w'-w) \\ &=& \iint\mathrm dx\mathrm dy\frac{\partial(x,y)}{\partial(z,w)}f_{Z,W}(z,w)\delta(z'-z)\delta(w'-w) \\ &=& \iint\mathrm dx\mathrm dyf_{X,Y}(x,y)\delta(z'-z(x,y))\delta(w'-w(x,y))\;. \end{eqnarray}