do you see anywhere in your book where this "DTFT" is defined for your w.s.s. process, $x[n]$? the DTFT is normally defined as
$$ X(e^{j \omega}) \triangleq \sum\limits_{n=-\infty}^{+\infty} x[n] e^{-j \omega n} $$
but for this infinite sum to converge, we normally require that
$$ \sum\limits_{n=-\infty}^{+\infty} \Big|x[n]\Big| < +\infty $$
which is a stronger requirement than
$$ \sum\limits_{n=-\infty}^{+\infty} \Big|x[n]\Big|^2 < +\infty $$
which says simply that the energy of $x[n]$ is finite, and $x[n]$ is a finite energy signal or just an energy signal.
a finite power signal or just a power signal does not satisfy the restrictions above, but satisfies this:
$$ \lim_{M \to +\infty} \frac{1}{2M+1} \sum\limits_{n=-M}^{+M} \Big|x[n]\Big|^2 < +\infty $$
sometimes, to do some math, we have to express this a little more strictly as
$$ \lim_{M \to +\infty} \frac{1}{2M+1} \sum\limits_{n=-M}^{+M} \Big|x[n]\Big| < +\infty $$
both restrictions are equivalently true, but the first one is expressing power directly.
so these "energy signals" and "power signals" live in different metric spaces, or more particularly, they live in different Hilbert spaces.
Hilbert spaces have something that we call "inner product" (or sometimes "dot product"):
$$ \big\langle x,y \big\rangle \quad \text{ or } \quad \big\langle x[n],y[n] \big\rangle $$
and we always define the norm of $x[n]$ in terms of this inner product:
$$ \big\| x \big\| \triangleq \sqrt{\big\langle x,x \big\rangle} $$
now the inner product for the space of energy signals is defined to be:
$$ \big\langle x[n],y[n] \big\rangle \triangleq \sum\limits_{n=-\infty}^{+\infty} x[n] \cdot \overline{y[n]} $$
where $\overline{y[n]}$ is the complex conjugate of $y[n]$. FYI, we see that the inner product is "conjugate commutative":
$$ \big\langle x,y \big\rangle = \overline{\big\langle y,x \big\rangle} $$
it should be easy to see that the norm $\big\| x \big\|$ is just the square root of the total energy:
$$\begin{align}
\big\| x \big\|^2 &= \big\langle x,x \big\rangle \\
&= \sum\limits_{n=-\infty}^{+\infty} x[n] \cdot \overline{x[n]} \\
&= \sum\limits_{n=-\infty}^{+\infty} \Big|x[n]\Big|^2 \quad < \infty \\
\end{align} $$
similarly, the inner product for the space of power signals is defined to be:
$$ \big\langle x[n],y[n] \big\rangle \triangleq \lim_{M \to +\infty} \frac{1}{2M+1} \sum\limits_{n=-M}^{+M} x[n] \cdot \overline{y[n]} $$
again, the inner product is conjugate commutative:
$$ \big\langle x,y \big\rangle = \overline{\big\langle y,x \big\rangle} $$
and we see that the norm $\big\| x \big\|$ is just the square root of the total power:
$$\begin{align}
\big\| x \big\|^2 &= \big\langle x,x \big\rangle \\
&= \lim_{M \to +\infty} \frac{1}{2M+1} \sum\limits_{n=-M}^{+M} x[n] \cdot \overline{x[n]} \\
&= \lim_{M \to +\infty} \frac{1}{2M+1} \sum\limits_{n=-M}^{+M} \Big|x[n]\Big|^2 \quad < \infty \\
\end{align} $$
now for both energy and power signals, we can define the cross-correlation between two signals as:
$$ R_{xy}[m] \triangleq \big\langle x[n+m], y[n] \big\rangle $$
for energy signals it's
$$\begin{align}
R_{xy}[m] & = \big\langle x[n+m], y[n] \big\rangle \\
& = \sum\limits_{n=-\infty}^{+\infty} x[n+m] \cdot \overline{y[n]} \\
\end{align}$$
for power signals it's
$$\begin{align}
R_{xy}[m] & = \big\langle x[n+m], y[n] \big\rangle \\
& = \lim_{M \to +\infty} \frac{1}{2M+1} \sum\limits_{n=-M}^{+M} x[n+m] \cdot \overline{y[n]} \\
\end{align}$$
the autocorrelation of $x[n]$, denoted $R_{xx}[m]$, is simply the cross-correlation of $x[n]$ with itself
$$ R_{xx}[m] = \big\langle x[n+m], x[n] \big\rangle $$
the spectral density is the DTFT of the autocorrelation
$$ S_{xx}(e^{j\omega}) = \sum\limits_{m=-\infty}^{+\infty} R_{xx}[m] \, e^{-j \omega m} $$
and, more generally, the cross-spectral density is the DTFT of the cross-correlation
$$ S_{xy}(e^{j\omega}) = \sum\limits_{m=-\infty}^{+\infty} R_{xy}[m] \, e^{-j \omega m} $$
now, for a deterministic (that means $x[n]$ is known for all $n$) and finite energy signal, the DTFT of $x[n]$ is defined
$$ X(e^{j \omega}) = \sum\limits_{n=-\infty}^{+\infty} x[n] e^{-j \omega n} $$
and it can be shown that the spectral density is simply the magnitude-square of the DTFT:
$$ S_{xx}(e^{j\omega}) = \Big| X(e^{j \omega}) \Big|^2 $$
the proof is as
$$\begin{align}
\Big| X(e^{j \omega}) \Big|^2 &= X(e^{j \omega}) \cdot \overline{X(e^{j \omega})} \\
&= \left( \sum\limits_{n=-\infty}^{+\infty} x[n] e^{-j \omega n} \right) \overline{\left( \sum\limits_{n=-\infty}^{+\infty} x[n] e^{-j \omega n} \right)} \\
&= \left( \sum\limits_{n=-\infty}^{+\infty} x[n] e^{-j \omega n} \right) \left( \sum\limits_{n=-\infty}^{+\infty} \overline{x[n] e^{-j \omega n}} \right) \\
&= \left( \sum\limits_{n=-\infty}^{+\infty} x[n] e^{-j \omega n} \right) \left( \sum\limits_{n=-\infty}^{+\infty} \overline{x[n]} e^{+j \omega n} \right) \\
&= \left( \sum\limits_{m=-\infty}^{+\infty} x[m] e^{-j \omega m} \right) \left( \sum\limits_{n=-\infty}^{+\infty} \overline{x[n]} e^{+j \omega n} \right) \\
&= \sum\limits_{m=-\infty}^{+\infty} \sum\limits_{n=-\infty}^{+\infty} x[m] e^{-j \omega m} \overline{x[n]} e^{+j \omega n} \\
&= \sum\limits_{m=-\infty}^{+\infty} \sum\limits_{n=-\infty}^{+\infty} x[m] \overline{x[n]} e^{-j \omega (m-n)} \\
&= \sum\limits_{m=-\infty}^{+\infty} \sum\limits_{n=-\infty}^{+\infty} x[n+m] \overline{x[n]} e^{-j \omega (n+m-n)} \\
&= \sum\limits_{m=-\infty}^{+\infty} \left(\sum\limits_{n=-\infty}^{+\infty} x[n+m] \overline{x[n]} \right) e^{-j \omega m} \\
&= \sum\limits_{m=-\infty}^{+\infty} R_{xx}[m] e^{-j \omega m} \\
&= S_{xx}(e^{j \omega}) \\
\end{align}$$
this is the case for a finite energy signal. but for a power signal, $x[n]$, the DTFT $X(e^{j \omega})$, does not simply exist without some hand-waving using the dirac delta function $\delta(\omega)$ in the frequency domain, which i do not want to do, and will stay away from that kind of representation. but the power spectral density $S_{xx}(e^{j \omega})$ continues to exist for power signals, even if the DTFT and $\Big|X(e^{j \omega})\Big|^2$ do not.
Filtering
whether these are both energy signals or both power signals, with $x[n]$ as input and $y[n]$ as output, a discrete-time LTI system (or "filter") will relate the two with the convolution summation:
$$\begin{align}
y[n] &= \sum\limits_{m=-\infty}^{+\infty} h[m] x[n-m] \\
&= \sum\limits_{m=-\infty}^{+\infty} h[n-m] x[m] \\
\end{align}$$
in the frequency domain, taking the DTFT of both sides is
$$ Y(e^{j \omega}) = H(e^{j \omega}) \cdot X(e^{j \omega}) $$
the impulse response $h[n]$ must be a finite energy signal
$$ \sum\limits_{n=-\infty}^{+\infty} \Big|h[n]\Big|^2 < +\infty $$
and the frequency response $H(e^{j \omega})$ is the DTFT of the impulse response
$$ H(e^{j \omega}) = \sum\limits_{n=-\infty}^{+\infty} h[n] e^{-j \omega n} $$
if $x[n]$ and $y[n]$ are energy signals, it follows directly that
$$\begin{align}
Y(e^{j \omega}) &= H(e^{j \omega}) \cdot X(e^{j \omega}) \\
\\
\Big|Y(e^{j \omega})\Big| &= \Big|H(e^{j \omega})\Big| \ \Big|X(e^{j \omega})\Big| \\
\\
\Big|Y(e^{j \omega})\Big|^2 &= \Big|H(e^{j \omega})\Big|^2 \ \Big|X(e^{j \omega})\Big|^2 \\
\\
S_{yy}(e^{j \omega}) &= \Big|H(e^{j \omega})\Big|^2 \ S_{xx}(e^{j \omega}) \\
\end{align}$$
but to show this for finite power signals is a little more difficult. remember that for $y[n]$, $S_{yy}(e^{j \omega})$ continues to exist, even if the DTFT $Y(e^{j \omega})$ and $\Big|Y(e^{j \omega})\Big|^2$ do not. and the convolution operation continues to exist with the same summation expression for power signals as with energy signals.
recall for power signals, the autocorrelation is:
$$\begin{align}
R_{yy}[m] & = \big\langle y[n+m], y[n] \big\rangle \\
& = \lim_{M \to +\infty} \frac{1}{2M+1} \sum\limits_{n=-M}^{+M} y[n+m] \cdot \overline{y[n]} \\
\end{align}$$