Here is my understanding of it. I will take the help of two widely read books: Griffiths QM - 2nd ed. (Ref. 1) and Sakurai's Modern Quantum Mechanics - revised ed. (Ref. 2)
PART 1: WHY ARE ELECTRONS CONSIDERED TO HAVE A SPIN ANGULAR MOMENTUM?
Since electrons respond to a magnetic field, they are said to have a magnetic moment, and a Hamiltonian $H = - \vec{\mu} \cdot \vec{B}$ is lifted from classical physics, with $\mu$ being the magnetic moment. Under this $H$, the e- have been found to rotate (or precess; i.e. Larmor precession). It is simply a precession of the magnetic moment which is detectable when placed between the Stern-Gerlach magnets. This precession has been worked out in Example 4.3 of Ref. 1
Sec. 3.2 of Ref. 2. Only if we assume $\vec{\mu} = \gamma \vec{S}$ ($\gamma$ = gyromagnetic ratio) the evolution operator $\exp(- i H t/ \hbar )$ becomes $\exp( i \vec{S} \cdot \vec{B} ~\gamma t/ \hbar )$ which is a rotation operator (recall that spin is the generator of rotations). Hence, ascribing the magnetic moment to some kind of angular momentum seems reasonable.
PART 2: HOW TO DEFINE/MEASURE THE MAGNITUDE OF SPIN?
We will work things out for the case of $s=1/2$, i.e. an e-. Higher $s$ cases can be dealt with similarly. Sec. 3.5 of Ref. 2 works out the eigenvalues of spin operators with
$$ [S_i,S_j] = \alpha ~\epsilon_{ijk} S_k, $$
as the starting point (the book takes $\alpha= i \hbar$ from the get-go but we will arrive at this conclusion later). These relations give the eigenvalues of $S^2$ and $S_z$ to be $-s(s+1) \alpha^2 $ and $-i m \alpha$ where $s$ is determined by the number $N$ into which a spin multiplet splits when the e- beam is sent through the Stern-Gerlach apparatus such that $N = 2 s+1$. All this comes from the angular momentum algebra presented in Sec. 3.5 of Ref. 2.
Further from Example 4.3 of Ref. 1, we can see that
$$
H = - \vec{\mu} \cdot \vec{B} , \\
\vec{\mu} = \gamma \vec{S} , \\
\omega = \gamma B ,
$$
where $\omega$ stands for the Larmor precession frequency. Upon measuring $\omega$ and $B$, $\gamma$ can be computed. Since the energy change associated with the change of one angular momentum/spin state to another can be measured, $\Delta H = - \Delta \vec{\mu} \cdot \vec{B} = - \gamma \Delta \vec{S} \cdot \vec{B} = - \gamma i (\Delta m) \alpha B$ lets us determine $\alpha$. We assume here that $\vec{B}$ points in the z-direction. My guess is that $\alpha$ comes out to be $i \hbar$ (which all textbooks assume from the get-go and write $ [S_i,S_j] = i \hbar~\epsilon_{ijk} S_k, $ instead of $ [S_i,S_j] = \alpha ~\epsilon_{ijk} S_k, $).
Hence, measured values (eigenvalues) of $S^2$ and $S_z$ become $s(s+1) \hbar^2$ and $m \hbar$, respectively. Note that to measure quantum spin, we did not measure some position and momentum vector and took its cross product. Associating the spin with a magnetic moment and studying its Larmor frequency as well as borrowing concepts from classical physics is crucial to operational define spin in quantum mechanics.
For more insight, check these out too:
Looking Glass Universe video 1
Looking Glass Universe video 2