Over the last couple of decades, the theory of interacting particle systems has found some unexpected connections to orthogonal polynomials, symmetric functions, and various combinatorial structures. The asymmetric simple exclusion process (ASEP) has played a central role in this connection. Recently, Cantini, de Gier, and Wheeler found that the partition function of the multispecies ASEP on a circle is a specialization of a Macdonald polynomial $P_{\lambda}(X;q,t)$. Macdonald polynomials are a family of symmetric functions that are ubiquitous in algebraic combinatorics and specialize to or generalize many other important special functions. Around the same time, Martin gave a recursive formulation expressing the stationary probabilities of the ASEP on a circle as sums over combinatorial objects known as multiline queues, which are a type of queueing system. Shortly after, with Corteel and Williams we generalized Martin's result to give a new formula for $P_{\lambda}$ via multiline queues.

The modified Macdonald polynomials $\widetilde{H}_{\lambda}(X;q,t)$ are a version of $P_{\lambda}$ with positive integer coefficients. A natural question was whether there exists a related statistical mechanics model for which some specialization of $\widetilde{H}_{\lambda}$ is equal to its partition function. With Ayyer and Martin, we answer this question in the affirmative with the multispecies totally asymmetric zero-range process (TAZRP), which is a specialization of a more general class of zero range particle processes. We introduce a new combinatorial object in the flavor of the multiline queues, which on one hand, expresses stationary probabilities of the mTAZRP, and on the other hand, gives a new formula for $\widetilde{H}_{\lambda}$. We define an enhanced Markov chain on these objects that lumps to the multispecies TAZRP, and then use this to prove several results about particle densities and correlations in the TAZRP.