Quasi-isometric rigidity of Nil
Here's a sketch of proof (which intersects yours):
Let $\Gamma$ be QI to NIL. As you say, by Gromov's theorem, $\Gamma$ is virtually nilpotent; let some finite index subgroup be a lattice in some simply connected nilpotent Lie group $G$.
The growth of $G$ and NIL are equivalent, hence the growth of $G$ is $\simeq r^4$. Only two simply connected nilpotent Lie groups have growth $\sim r^4$: the abelian 4-dimensional one, and NIL. (This is a baby case of results of Guivarch or Jenkins in 1973)
So either $\Gamma$ is virtually $\mathbf{Z}^4$, or it's virtually a lattice in NIL.
To exclude the 1st case, you can invoke Pansu's 1983 paper but it's indeed a baby case: it's enough to show that NIL and $\mathbf{R}^4$ are not QI. But, NIL being endowed with a Carnot-Catheodory metric, both spaces are proper spaces with nontrivial 1-parameter groups of non-isometric similarities, so are isometric to their asymptotic cones. Since they are not homeomorphic, it follows that the asymptotic cones are not isometric and hence NIL and $\mathbf{R}^4$ are not quasi-isometric. Avoiding asymptotic cones, another approach is the Dehn function, which is $\simeq n^2$ for $\mathbf{R}^4$ and $\simeq n^3$ for NIL. Also, Roe proved later that for a simply connected solvable Lie group, the dimension is a quasi-isometry invariant.
Note: While's it's not seriously needed here, Pansu's 1983 result is in general much harder because it describes the asymptotic cone for simply connected nilpotent Lie groups that are not Carnot, i.e., do not admit a Carnot-Caratheodory metric with a non-isometric 1-parameter group of dilations. It gives a description of the asymptotic cone up to bilipschitz homeomorphism, but to get strong consequences one has to combine with Pansu's 1989 Annals result, which says that the associated Carnot-graded Lie algebra (a purely algebraic object! no metric here) of a simply connected nilpotent Lie group is a QI-invariant. This allows to distinguish, for instance, some infinite families in fixed dimension.
Edit: there is a statement avoiding passing to finite index subgroups: if $\Gamma$ is a finitely generated group QI to NIL, (NIL is being endowed with a left-invariant metric invariant under a maximal compact subgroup of automorphisms), then $\Gamma$ admits a proper cocompact isometric action on NIL.
Here's another, more topological way to see that $\mathrm{Nil}$ is not quasi-isometric to $\mathbf{R}^4$: use the fact that $\mathrm{Nil}$ is homeomorphic to $\mathbf{R}^3$, so that $\pi_2^\infty(\mathrm{Nil})$ is nontrivial, whereas $\pi_2^\infty(\mathbf{R}^4)$ is trivial.
This works because $\pi_2^\infty$ is indeed a quasi-isometry invariant within a class of metric spaces which contains $\mathbf{R}^4$ and $\mathrm{Nil}$. The key property is that they are uniformly 2-connected, i.e. spheres of dimension $\le 2$ can be filled by balls with metric control.
Concretely, to prove the result by hand you can do the following: assume that there is a quasi-isometry $f:\mathrm{Nil} \to \mathbf{R}^4$. Without loss of generality, $f$ takes the origin to the origin. Let $\phi$ be an embedding of $\mathbf{S}^2$ into $\mathrm{Nil}$ with image a metric sphere of large radius centred at the origin. Put a sufficiently fine triangulation $\mathcal T$ on $\mathbf{S}^2$. By induction on the skeleta, construct a continuous map $g:\mathbf{S}^2\to \mathbf{R}^4$ which coincides with $f\circ\phi$ on the $0$-skeleton and is at bounded distance from $f\circ\phi$.
Since $\mathbf{R}^4$ is 2-connected at infinity, we can extend $g$ to a continuous map $\bar g$ from the 3-ball to $\mathbf{R}^4$ which stays far away from the origin. Composing with a coarse inverse of $f$ and using the same trick to make it continuous, we obtain a continous map from the $3$-ball to $\mathrm{Nil}\setminus\{O\}$ which extends $\phi$. This is a contradiction.
(Remark: we have actually only used that $\mathrm{Nil}$ is uniformly 2-connected and $\mathbf{R}^4$ is uniformly 1-connected.)