This answer nicely categorises which permutations of a finite set have square roots. This prompts the following question:
Does every permutation $\sigma:\mathbb{N}\to\mathbb{N}$ that contains no finite cycles (i.e. it is "finite-cycle-free") have a square root?
Finite-cycle-free permutations can be generated by using a double (two-way) greedy algorithm that avoids cycles.
For example, $\sigma(1) = 2,\ \sigma^{-1}(1) = 3,\ \sigma(2) = 4,\ \sigma^{-1} (3) = 5,\ \sigma(4) = 6,\ \sigma^{-1} (5) = 7,\ldots. $
You can generate other finite-cycle-free permutations by choosing larger numbers, in other words, the algorithm doesn't have to be greedy. For example, $\sigma(1) = 12,\ \sigma^{-1}(1) = 27,\ \sigma(2) = 400,\ \ldots. $ But you must avoid cycles.
Anyway, the question above is what I am interested in. I haven't found a square root for the greedy algorithm in the first example. So maybe there isn't one, but right now I don't see a reason why this would be the case.
Any cycle-free permutation of $\mathbb{N}$ decomposes $\mathbb{N}$ into a bunch of infinite chains, each of which "looks like" $\mathbb{Z}$ equipped with the successor map; whether or not a square root exists is appropriately isomorphism-invariant, so the only relevant variable is how many orbits the permutation has.
If $\sigma$ has an even number of orbits, then $\sigma$ has a square root. This is because we can just "bounce between" pairs of orbits. For instance, if $\sigma$ has two orbits, then - reconstruing $\sigma$ as acting on $\mathbb{Z}\times\{0,1\}$ via $\sigma(z,i)=(z+1, i)$ - we can get a square root via $$(z,0)\mapsto (z,1), \quad (z,1)\mapsto (z+1,0).$$ And this generalizes easily to any even number of orbits.
It remains to non-square-rootability of the permutations with an odd number of orbits. I strongly suspect this is true, but I don't immediately see a "clean" proof of this. For example, let's consider the case of a single orbit, for simplicity thinking of $\mathbb{Z}$ equipped with the successor map $\sigma: z\mapsto z+1$. Assume $\tau^2=\sigma$. We have for example $$\tau(0)=2\implies \tau(2)=1\implies \tau(1)=3\implies\tau(3)=2\implies\tau(2)=4,$$ so $\tau(0)$ can't be $2$. Similar arguments rule out $\tau(0)$ being any number whatsoever, so indeed the one-orbit cycle-free permutation on $\mathbb{N}$ has no square root. But generalizing this to arbitrarily many odd orbits seems messy.