A recent paper \cite{KMMO} introduced the stochastic U_q(A_n^{(1)}) vertex model. The stochastic S-matrix is related to the R-matrix of the quantum group U_q(A_n^{(1)}) by a gauge transformation. We will show that a certain function D^+_μ intertwines with the transfer matrix and its space reversal. When interpreting the transfer matrix as the transition matrix of a discrete-time totally asymmetric particle system on the one-dimensional lattice Z, the function D^+_μ becomes a Markov duality function D_μ which only depends on q and the vertical spin parameters μ_x. By considering degenerations in the spectral parameter, the duality results also hold on a finite lattice with closed boundary conditions, and for a continuous-time degeneration. This duality function had previously appeared in a multi-species ASEP(q,j) process. The proof here uses that the R-matrix intertwines with the co-product, but does not explicitly use the Yang-Baxter equation. It will also be shown that the stochastic U_q(A_n^{(1)}) is a multi-species version of a stochastic vertex model studied in \cite{BP,CP}. This will be done by generalizing the fusion process of \cite{CP} and showing that it matches the fusion of \cite{KRL} up to the gauge transformation. We also show, by direct computation, that the multi-species q-Hahn Boson process (which arises at a special value of the spectral parameter) also satisfies duality with respect to D_0, generalizing the single-species result of \cite{C}.