Spin operator in an arbitrary direction let $\ \hat{n} = (\sin{\phi}\cos{\theta},\sin{\phi} \sin{\theta},\cos{\phi})$
How to prove $\vec{S}\cdot\hat{n}$ is the spin operator in an arbitrary axis rigorously?
 A: You just have to find an eigenvalues of this operator. Suppose you have the case of spin $\frac{1}{2}$ particle. Then
$$
\vec{S} = \frac{1}{2}\left( \sigma_{x}, \sigma_{y}, \sigma_{z}\right),
$$
where $\sigma_{i}$ is Pauli matrices.
Now you have that 
$$
\hat{h}\equiv \vec{S}\cdot \vec{n} =\sum_{\sigma = \pm \frac{1}{2}}\sigma| v_{\sigma}\rangle\langle v_{\sigma} | ,
$$
where $v_{\sigma}$ is the eigenvector of $\hat{h}$ operator (which, by the correspondence principle, defines projection of spin on $\vec{n}$ axis) and $\sigma$ is helicity. Physical sense of eigenvector is that expansion coefficient near it for an arbitrary state $|\psi \rangle$,
$$
|\psi\rangle = \sum_{\sigma}c_{\sigma}|v_{\sigma}\rangle, \quad \sum_{\sigma}|c_{\sigma}|^{2} = 1
$$ 
defines the probability to obtain value $\sigma$ when we measure average value of $\hat{h}$,
$$
\langle \psi|\hat{h} | \psi\rangle =\sum_{\sigma}\sigma |c_{\sigma}|^{2}
$$
It is clear that $|v_{\sigma}\rangle$ is the the state in which the value of the spin in the direction given by $\vec{n}$. Really, by setting $\theta = \varphi = 0$ you obtain that $\hat{h} = \sigma_{z}$, while $v_{\sigma}= \begin{pmatrix} 1 \\ 0\end{pmatrix}, \begin{pmatrix} 0\\ 1\end{pmatrix}$ are eigenvectors of $\sigma_{z}$.
