Justification for "reducing integrals" in the virial expansion for gases

I have been following this document (https://sites.chem.utoronto.ca/chemistry/jmschofi/chm427/gases.pdf) regarding the virial expansion of gases and on finding the virial coefficients.

On Page 7, they do the following:
$$Z_3 = int _V dmathbf{r}_1 dmathbf{r}_2 dmathbf{r}_3 f_{12}f_{23}$$
where $f_{ij} = f(|mathbf{r}_i-mathbf{r}_j|)$, the Meyer-f function. We are integrating over the entire phase space, with $mathbf{r}$ being the 3D position vector.
They say that the above integral is reducible, and reduce it as,
$$Z_3 = int _V dmathbf{r}_1 dmathbf{r}_2 dmathbf{r}_3 f_{12}f_{23} = int _V dmathbf{r}_{12} dmathbf{r}_{2} dmathbf{r}_{23} f_{12}f_{23}=Vleft(int _V dmathbf{r}_{12}f_{12}right)left( int _Vdmathbf{r}_{23}f_{23}right)$$

My question is, why is this last step valid?

I understand substitution (like $u$-substitution, or trig substitution when we are learning basic integrals) of variables within an integral sign, but this seems to be a bit more than that, right?
I made the substitution, $$mathbf{r}_1 rightarrow mathbf{r}_{12} mathbf{r}_{3} rightarrow mathbf{r}_{23}$$
But, if $r_{12} = r_1 – r_2$, then shouldn’t $dmathbf{r}_{12} = dmathbf{r}_1 – dmathbf{r}_2$? Furthermore, if $f_{12}$ and $f_{23}$ have a dependence on $mathbf{r}_2$, why are we ignoring that dependence and just integrating it out?

I am really struggling to understand why that step is allowed.

I would appreciate any advice you have for me.