Given a quaternion $q$ and a complex number $\lambda$, for scalar multiplication you can either apply $\lambda$ to $q$ from the left (i.e. $\lambda q$) or from the right (i.e. $q\lambda$). This means you can interpret $\mathbb{H}$ as a left complex vector space or as a right complex vector space. Now consider the left and right multiplication maps
$$ L_p(x)=px, \qquad R_p(x)=xp. \tag{1}$$
Then $L_p(x\lambda)=L_p(x)\lambda$ for all complex numbers, so $L_p$ is a linear transformation of $\mathbb{H}$ as a right complex vector space. And $R_p(\lambda x)=\lambda R_p(x)$, so $R_p$ is a linear transformation of $\mathbb{H}$ as a left complex vector space.
However, $L_p$ is not linear if we treat $\mathbb{H}$ as a left complex vector space, and $R_p$ is not linear if we treat $\mathbb{H}$ as a right vector space. This is because $\mathbb{H}$ is not commutative. (Exercise.)
Using $\{1,\mathbf{j}\}$ as a basis for $\mathbb{H}$ as a complex vector space, you would write an arbitrary quaternion as $z+w\mathbf{j}$ if you're thinking left vector space, and write $z+\mathbf{j}w$ if you're thinking right vector space. (Also notice $w$ comes before $z$ in the alphabet, so we're backwards alphabetically.) You're thinking left vector space but you're examining the right-linear transformation $L_p$. That's problematic.
If you look at $L_p$ and think right vector space, you should get
$$ \begin{bmatrix} z & -\overline{w} \\ w & \phantom{-}\overline{z} \end{bmatrix}. \tag{2} $$
(Exercise.)
Representing quaternions as matrices using $R_p$ is problematic since $R_p\circ R_q\ne R_{pq}$. Indeed, this would be a function $\mathbb{H}\to M_2(\mathbb{C})$ which is not a homomorphism but rather an anti-homomorphism, i.e. it satisfies $R_p\circ R_q=R_{qp}$ (the order of multiplication is reversed). In order to turn it into a homomorphism proper, one must either (pre)compose it with an anti-automorphism of $\mathbb{H}$ or (post)compose with an anti-automorphism of $M_2(\mathbb{C})$. Quaternion conjugation satisfies $\overline{pq}=\overline{q}\,\overline{p}$ and matrix transpose satisfies $(AB)^T=B^TA^T$, so we can use these as anti-automorphisms.
In the case of quaternion conjugation, we can take $p=z+w\mathbf{j}$, then its conjugate $\overline{p}=\overline{z}-w\mathbf{j}$, then the matrix of $R_{\overline{p}}$ you can calculate (exercise) to be
$$ \begin{bmatrix} \overline{z} & \overline{w} \\ -w & z \end{bmatrix}, \tag{3} $$
and if instead you just calculated $R_p$'s matrix and took the transpose you'd get $(2)$ again.
Complex conjugation (applied entry-wise to a matrix) is an automorphism of $M_2(\mathbb{C})$, i.e. it satisfies $\overline{AB}=\overline{A}\,\overline{B}$, so we may compose it with any of the representations above to get another valid representation.
Let's skip characteristic $2$ for the moment, so I'll assume $F$ is a field with characteristic $\ne2$.
The $F$-algebra $H(F)$ of quaternions over $F$ that you define has the antiautomorphism $q\mapsto q^*$, where, for $q=w+xi+yj+zk$, $q^*=w-xi-yj-yk$. It is bijective and preserves addition, but satisfies
$$
(q_1q_2)^*=q_2^*q_1^*
$$
by direct computation.
If we define $N(q)=qq^*=q^*q$, we see that $N(q)\in F$. Now, if $N(q)\ne0$, the quaternion $q$ is invertible, because
$$
\frac{q^*}{N(q)}q=1=q\frac{q^*}{N(q)}
$$
On the other hand, if $N(q)=0$, the quaternion $q$ is a zero-divisor by definition.
Suppose $H(F)$ is not a division ring. Take a nonzero noninvertible quaternion $q=w+xi+yj+zk$; then $N(q)=w^2+x^2+y^2+z^2=0$ and therefore we can write $-1$ as a sum of squares in $F$.
According to Artin-Schreier theory of formally real fields, a field where $-1$ is not a sum of squares can be ordered, so it has characteristic $0$ and is infinite.
Over a field of odd prime characteristic, it may happen that $-1$ is a sum of squares, but you need more than three. However, the four-square theorem you mention proves that the quaternion algebra is not a division ring also in this case: indeed, if the characteristic is $p$, then $p=a^2+b^2+c^2+d^2$ for some integers $a,b,c,d$ and therefore the quaternion $q=a+bi+cj+dk$ (where $a,b,c,d$ are interpreted in $F$, as usual) has $N(q)=0$.
By the way, in Herstein's “Topics in Algebra”, Wedderburn's theorem is exploited to prove the four-square theorem.
What happens in characteristic $2$? Well, the conjugation is the identity, so it's not really useful. However, the quaternion $1+i$ is not invertible, because $(1+i)^2=0$.
In conclusion, $H(F)$ is never a division ring, unless $F$ has characteristic $0$ and $-1$ cannot be written as a sum of at most three squares. Finiteness of $F$ is not needed.
Best Answer
The field of complex numbers is isomorphic to a subfield of the division ring of the quaternions. Also, the field of real numbers is isomorphic to a subfield of the complex field.
Whether or not we have subsets here depends on the way the complex numbers and the quaternions are defined.