It is well known that two randomly chosen permutations of $n$ symbols commute with probability $p_n/n!$ where $p_n$ is the number of partitions of $n$. This is a special case of the fact that in a group, the probability that two elements chosen uniformly at random (with repetition allowed) commute is the number of conjugacy classes divided by the size of the group.
Question. What is the probability that two mappings of $n$ symbols chosen uniformly at random commute?
I suspect an exact answer would be difficult and would be happy to learn of reasonably tight asymptotic results.
Added. This probability should go to zero quickly because Misha Berlinkov recently showed that with probability going to 1 as $n$ goes to infinity, two random elements generate a subsemigroup containing a constant map and so if they commute they generate a unique constant map. This should happen almost never (and most likely has been proven).
Added based on Brendan McKay's answer. Computing the probability that an element of a monoid $M$ commutes with an element of its groups of units $G$ is no harder than the commuting probability in a group. Namely, $G$ acts on $M$ by conjugation; let's call the orbits conjugacy classes. Then the probability that an element of $G$ commutes with an element of $M$ is the number of conjugacy classes of $M$ divided by the number of elements of $M$. The proof is the same as for groups. If $Fix(g)=\{m\in M\mid gmg^{-1}=m\}$, then
$$\frac{|\{(g,m)\in G\times M\mid gm=mg\}|}{|G||M|} = \frac{1}{|M|}\frac{1}{|G|}\sum_{g\in G}|Fix(g)| = \frac{\text{number of conjugacy classes}}{|M|}$$ by the Cauchy-Burnside-Frobenius orbit counting formula.
For $M=T_n$ the monoid of all mappings on $n$ symbols and $G=S_n$ the symmetric group, conjugacy classes correspond to isomorphism classes of functional digraphs on $n$ vertices. A functional digraph is a digraph (loops allowed) in which each vertex has outdegree $1$. Each mapping $f$ gives a functional digraph by drawing an edge from $i$ to $f(i)$. It is obvious that $f,g$ are conjugate iff their corresponding digraphs are isomorphic (it is the same proof that permutations are determined up to conjugacy by cycle type).
According to the book of Flajolet and Sedgewick, the number of unlabelled functional digraphs grows likes $O(\rho^{−n}n^{−1/2})$ where $\rho\approx .29224$. So the probability of a random mapping commuting with a random permutation is pretty small. Brendan raises the nice question of how different the probability of a random permutation commuting with a random mapping is from the probability of a random mapping commuting with a random mapping. My guess is the latter goes to $0$ qualitatively faster.
Best Answer
The number of ordered pairs of commuting functions is A181162. I agree with those counts up to n=7. There is little in OEIS that helps to answer the asymptotics question.
Incidentally, the probability that $f(g(1))=g(f(1))$ is not $1/n$. I think it is $1/n + (n-1)/n^3~$ though I might have miscalculated. That formula works up to $n=7$.
ANOTHER relevant fact: If $f$ is a permutation, then any function $g$ commuting with $f$ is determined by the image of one element of each cycle of $f$. So the number of such $g$ is at most $n^{C(f)}$ where $C(f)$ is the number of cycles of $f$. Random permutations have on average only $\ln n+O(1)$ cycles, so the probability of a random function commuting with a random permutation might be at most something like $n^{-n+\ln n+O(1)}$ (which is an abuse of expectations but might be something akin to the truth). Is a random function more or less likely to commute with another random function or with a random permutation? [NOTE: I added "at most" since some assignments don't work: the image of a point in a cycle of length $k$ must lie in a cycle whose length is a divisor of $k$.]