Dirac’s goal was to find a first order time derivative which worked with $E^2=c^2 |\vec p|^2+m^2 c^4$
For this he proposed a free particle Hamiltonian of the form
$$ \boxed{H=c\vec \alpha \cdot \vec p + \beta m c^2} $$
where $m$ is the rest mass, $\vec \alpha$ has 3 components, $\beta$ has one and $\vec \alpha$ and $\beta$ don't necessarily commute with each other but do commute with $\vec p$
Squaring the equation gives
$$ H^2 =c^2 (\vec \alpha \cdot \vec p )(\vec \alpha \cdot \vec p ) + mc^3 (\vec \alpha \cdot \vec p ) \beta +mc^3 \beta(\vec \alpha \cdot \vec p ) + m^2 c^4 \beta ^2 $$
Comparing with $E^2 = c^2 \vec p ^2 +m^2 c^4$ we have
$$ \begin{aligned} (\vec \alpha\cdot\vec p ) (\vec \alpha \cdot \vec p ) & = |\vec p| ^2 \\ (\vec \alpha \cdot \vec p )\beta + \beta (\vec \alpha \cdot \vec p ) = (\vec \alpha \beta +\beta\vec \alpha ) \cdot \vec p &=0 \\ \beta ^2 &=1
\end{aligned} $$
Which rewriting in terms of anticommutators as
$$ \{\alpha _i ,\alpha _j \} =2\delta _{ij} \qquad \{ \alpha_i ,\beta\} =0 \qquad \beta ^2 =1 $$
where $i,j=1,2,3$ and $\{a,b\} =ab+ba$
💎 Conclusion: $\alpha _i , \beta$ cant be numbers but they can be matrices
🧽 Assume: $\alpha_i$ and $\beta$ are Hermitian matrix so that the Hamiltonian is Hermitian
🗒️ Note: since $\beta ^2 =1 = \alpha _{i}^2$ they are unitary and have eigenvalue $\pm 1$
💼 Case: lets consider the following equation to find the dimensions of the matrices
$$ \alpha _i \beta + \beta \alpha _i = 0 $$
Multiplying from the right by $\beta$ we get
$$ \alpha _i = -\beta \alpha _i \beta $$
where we used $\beta^2=1$
Taking the trace we can show that
$$ \operatorname{tr} (\alpha_i) =-\operatorname{tr} (\beta \alpha _i \beta ) =-\operatorname{tr} (\beta ^2 \alpha_i) =-\operatorname{tr} (\alpha_i) $$
💎 Conclusion: thus $\operatorname{tr} (\alpha_i)=0$ and in a similar argument $\operatorname{tr} (\beta)=0$
Since the eigenvalues of $\alpha _i$ and $\beta$ are $\pm 1$ and that the trace of $\alpha_i ,\beta$ is zero then they have even dimensions ie ($n\times n )$ where $n$ is even
Lets now define it as 4 dimensional as follows where we use the Dirac representation
$$ \vec \alpha =\begin{bmatrix}
0 & \vec \sigma \\ \vec \sigma & 0 \end{bmatrix} \qquad \beta = \begin{bmatrix} I_2 & 0 \\ 0 & -I_2 \end{bmatrix} $$
We can now replace $H$ and $\vec p$ by their operators and add a wavefunction $\psi ( \vec r,t)$ to give
$$ i\hbar \frac{\partial}{\partial t} \psi ( \vec r,t)=(c\vec \alpha \cdot \hat{\vec p} +\beta m c^2 ) \psi ( \vec r,t )=(-i\hbar c\vec \alpha \cdot \vec \nabla + \beta mc^2 ) \psi ( \vec r,t) $$
🗒️ Note: here the dot product is
$$ \vec \alpha \cdot \hat{\vec p } =\begin{bmatrix} 0 & \vec \sigma \cdot \hat{\vec p} \\ \vec \sigma \cdot \hat{\vec p} & 0 \end{bmatrix} \qquad \text{where}\qquad \vec \sigma \cdot \hat{\vec p} =\begin{bmatrix} \hat p_z & \hat p_x-i\hat p_y \\ \hat p_x + i\hat p_y & -\hat p_z \end{bmatrix} $$
Because the matrices are all four dimensional then so is the wave function
$$ \psi ( \vec r,t) = \begin{pmatrix} \psi_1 ( \vec r,t) \\ \psi_2 ( \vec r,t) \\ \psi_3 ( \vec r,t) \\ \psi_4 ( \vec r,t) \end{pmatrix} $$
🗒️ Note: this is not a relativistic four-vector instead it is a four-spinor
💎 Conclusion: the Dirac equation allows relativistic description of particles of spin-$\frac 12$