The Hubbard Model is already manifestly SU(2) spin-rotation invariant because there are no terms that connect the $\uparrow$ and $\downarrow$ sectors (here by connect, I mean change a $\uparrow$ to a $\downarrow$ or vice versa. The full Hamiltonian in terms of fermion annihilation/creation operators is:
$$ H = - t \sum \limits_{t,\sigma} \left( c^\dagger_{i,\sigma} c_{i+1,\sigma} + c^\dagger_{i+1,\sigma} c_{i,\sigma} \right) + U\sum_i n_{i,\uparrow} n_{i,\downarrow}\tag{1}$$
Notice that the hopping term only acts on either the $\uparrow$ or $\downarrow$ sectors, it never changes the **number$$ of $\uparrow$ or $\downarrow$ electrons. The potential ($U$) term does allow interaction between $\uparrow$ and $\downarrow$ electrons, but it is diagonal, so it never changes the total spin of the system.
In terms of spin ladder operators:
The most common way of thinking about the Hubbard model in terms of spin ladder operators is to consider the case of large $U$ at half filling. Half filling means there are $N$ electrons for $N$ sites, and large-$U$ means that the energy cost of double-occupying a site is high.
From this setup, we can do time-independent degenerate perturbation theory, treating the $U$ term as $H_0$ and the hopping as a perturbation. Let's consider ju From this standpoint, we have a set of degenerate ground states with one electron on each site. The energy of this ground state is zero ($E_0 =0$) and all possible combinations of up and down spins are degenerate.
As a rough example, let's consider an isolated pair of Hubbard sites with hamiltonian $H_0$ and hopping perturbation $V$:
\begin{eqnarray}
H_0 &= & U ( n_{+,1} n_{-,1} + n_{+,2} n_{-,2} ) \tag{2}\\
V &= & -t ( c^\dagger_{\pm,1} c_{\pm,2} + c^\dagger_{\pm,2} c_{\pm,1} )\tag{3}
\end{eqnarray}
The eigenstates of $H_0$ in the subspace of half filling are
\begin{eqnarray}
H_0 | 1 \rangle = & H_0 | \pm, \pm \rangle = 0 \tag{4}\\
H_0 |2 \rangle = & H_0 |\pm, \mp \rangle = 0 \tag{5}\\
H_0 |3 \rangle = & H_0 |0,2\rangle = U \tag{6}\\
H_0 |4 \rangle = & H_0 |2,0\rangle= U\tag{7}
\end{eqnarray}
This hamiltonian has two degenerate ground states: $|1\rangle,|2\rangle$ with eigenvalues of 0. Here, first order degenerate perturbation theory fails because there the matrix elements $\langle 1|V|2 \rangle = 0$, so we must advance to second-order perturbation theory:
$$
E^{(2)}_i = \sum \limits_{i \ne j} \frac{ \langle i | H_0 | j \rangle \langle j | H_0 | i \rangle }{E^{(0)}_i - E^{(0)}_j}\tag{8}
$$
In second order, a pair of antiparallel spins can make a virtual transition to a double occupied state and back. This effectively lowers the energy of antiparallel neighbors, while parallel neighbors don't benefit from this virtual transition (also called super exchange). For example, consider state 2:
\begin{align}
E^{(2)}_2 &= \frac{ \langle \pm \mp | H_0 | 0 2 \rangle \langle 0 2 | H_0 | \pm \mp \rangle }{0 - U} + \frac{ \langle \pm \mp | H_0 | 2 0 \rangle \langle 2 0 | H_0 | \pm \mp \rangle }{0 - U} \tag{9}\\ \nonumber
E^{(2)}_2 &= -\frac{t_0^2 }{U} - \frac{ t_0^2 }{ U} \tag{10}\\
&= - \frac{ 2 t_0^2}{U} \tag{11}
\end{align}
This means that $H |2 \rangle = -\frac{2t_0^2}{U} | 2\rangle $ therefore $J = \frac{2t_0^2}{U}$. This is not a rigorous derivation, but it gets the basic idea across.
The equivalent spin-ladder operator Hamiltonian to the large-U half-filling Hubbard model is there for just the Heisenberg model with $J = 2\frac{t^2}{U}$.
