@MattL 's link solves your problem, but it seems you need more help. Then I would like to provide here the same solution in different approach, expecting it might help.
Analysis of a downsampler block (more specifically the decimator or the compressor block) can be done in a number of ways. I prefer the following approach which I find it easier.
The first step of the analysis begins by observing that the structure of the figure a and figure b produce exactly the same outputs $y[n]$ given the same inputs $x[n]$, hence they are equivalent:
$$ \boxed{ x[n] \rightarrow \boxed{\downarrow M} \rightarrow y[n] ~~~~~~,~~~~~~~~~~~ \text{figure-a} }$$
$$ \boxed{
x[n] \rightarrow \boxed{ (\times)} \rightarrow v[n] \rightarrow \boxed{\downarrow M} \rightarrow y[n] ~,~ \text{figure-b} \\
~~~~~~~~~~~~~~~\uparrow w[n] }$$
Where the multiplier $w[n]$ is a periodic impulse train with a period of $M$ samples, which constitudes a discrete-time counterpart of a continuous-time impulse train sampler.
$$ w[n] = \sum_k \delta [n- kM] $$
Such kind of a structural modification, enables an unambiguous mathematical analysis at the expense of an increased number of steps and structural complexity. We shall analyse the figure-b and map its result to figure-a which is what we want indeed.
Then, as the second step, express the periodic impulse train $w[n]$, in the figure-b, by its (DFS) discrete-Fourier series sum:
$$ w[n] = \frac{1}{M} \sum_{k=0}^{M-1} e^{j \frac{ 2\pi k}{M}n }$$
If you cannot see why this is so, consult a DSP textbook for DFS.
In the third step observe that you can treat $v[n]$ as an expanded version of $y[n]$. Based on this, express Z-transforms of the expanded signal $v[n]$ in terms of Z-transform of $y[n]$ which is a lot easier to derive than the inverse operation of downsampling; more specifically:
$$V(z) = Y(z^M) \longrightarrow Y(z) = V(z^{\frac{1}{M}})$$
In the fourth step, we shall express Z-transform $V(z)$ of $v[n]$ in terms of Z-transform $X(z)$ of $x[n]$, based on the following two properties of Z-transforms:
$$ \sum_k x_k[n] \longrightarrow \sum_k X_k(z) ~~~ , ~~~ \text{linearity}$$
$$ e^{j\omega_0 n} x[n] \longrightarrow X(e^{-j \omega_0} z) ~~~ , ~~~ \text{complex modulation}$$
Then express $v[n] = x[n] w[n]$ and apply the above two properties:
$$ v[n] = x[n]w[n] = x[n] \frac{1}{M} \sum_{k=0}^{M-1} e^{j \frac{ 2\pi k}{M}n } = \frac{1}{M} \sum_{k=0}^{M-1} e^{j \frac{ 2\pi k}{M}n } x[n] $$
$$ V(z) = \frac{1}{M} \sum_{k=0}^{M-1} X( e^{-j \frac{ 2\pi k}{M}} z) $$
In the final step, merge the above result with the result of step-3 by replacing $z$ with $e^{j\omega}$ to convert the result from Z-transform to discrete-time Fourier transform:
$$ Y(z) = V(z^{1/M})$$
$$ V(z) = \frac{1}{M} \sum_{k=0}^{M-1} X( e^{-j \frac{ 2\pi k}{M}} z) $$
$$ Y(z) = \frac{1}{M} \sum_{k=0}^{M-1} X( e^{-j \frac{ 2\pi k}{M}} z^{1/M}) $$
$$ Y(z) = \frac{1}{M} \sum_{k=0}^{M-1} X( e^{-j \frac{ 2\pi k}{M}} e^{j \omega/M}) $$
which is the relation between the discrete time Fourier transforms of the input $x[n]$ of a compressor and its output $y[n] = x[Mn]$ as:
$$ Y(e^{j\omega}) = \frac{1}{M} \sum_{k=0}^{M-1} X( e^{j \frac{w- 2\pi k}{M}}) $$
Interpretation of the final result is as follows: The discrete time Fourier transform $Y(e^{j\omega})$ of the the downsampled signal $y[n] = x[Mn]$ is given by the shifted and scaled sum of $M$ copies of DTFT $X(e^{j\omega})$ of $x[n]$ scaled in magnitude by $1/M$ and expanded in frequency by $M$ and each copy shifted by $2\pi$ to the right.
This expansion in frequency $\omega$ by M of the DTFT $X(e^{j\omega})$ can be pictured as follows: Consider the original input frequency spectrum $X(e^{j\omega})$, which is bandlimited to $\omega = \pi/M$ for preventing aliasing due to subsequent downsampling. Being a DTFT, by default, It's periodic by $2\pi$ having its centers at $\omega = ... -2\pi, 0 , 2\pi, 4\pi, 6\pi ...$. Now expanding the frequency axis by M yields a function $X(e^{j\omega/M})$ such that it's periodic by $2 M \pi$ and it has its centers at $\omega = ..., -2M\pi, 0 , 2M\pi, 4M\pi, 6M\pi ...$. This spectrum has complete zeros at those frquencies $\omega = 2\pi, 4\pi, (M-1)2\pi ...$, due to frequency axis scaling, and those M copies of these expanded spectrums are added to each other after being shifted by $2\pi k$ to the right, for $k=0,1,2,M-1$, and therefore filling the gaps in between $\omega = 0$ and $\omega = M2\pi$ of the spectrum of $X(e^{j\omega/M})$, hence having a period of $2\pi$, $Y(e^{j\omega})$ is now a proper DTFT.
