I'm not sure about what "functorial" would entail here, but at least
when $p\ge5$ from a naive point of view things are quite simple. Once
one knows that the $j$-invariant of a supersingular
elliptic curve $E$ lies in $k=\mathbb{F}_{p^2}$ then there is a curve $E'$ 
defined over $k$ with the same $j$-invariant as $E$.
Up to $k$-isomorphism there are two candidates
for $E'$ but they are quadratic twists: one has $(p+1)^2$ points
over $k$ and the other has $(p-1)^2$ points. Equivalently
the $k$-Frobenius acts on one as $-p$ and the other as $+p$.

Let's pick an isomorphism $\alpha:E\to E'$ for each supersingular
curve where $E'$ is defined over $k$ with Frobenius $-p$. Given an
isogeny $\phi:E_1\to E_2$ then there is a corresponding
isogeny $\phi':E_1'\to E_2'$ making the obvious square commute.
This isogeny $\phi'$ is defined over $k$, because
it commutes with the $k$-Frobenius which on both sides equal $-p$.

We could proceed in exactly the same way taking each $E'$ to have
Frobenus $+p$ and come to the same conclusion. In some sense though, 
choosing $-p$ is more natural. If $E$ has $j$-invariant in
$\mathbb{F}_p$ then $E'$ will be defined over $\mathbb{F}_p$
if we take the $-p$ option but not the $+p$ option.

Of course in characteristic $2$ and $3$ things are different,
but in each case there is only one supersingular $j$-invariant.