Two-sided permutation statistics via symmetric functions

Given a permutation statistic $\operatorname{st}$, define its inverse statistic $\operatorname{ist}$ by $\operatorname{ist}(\pi):=\operatorname{st}(\pi^{-1})$. We give a general approach, based on the theory of symmetric functions, for finding the joint distribution of $\operatorname{st}_{1}$ and $\operatorname{st}_{2}$ whenever $\operatorname{st}_{1}$ and $\operatorname{st}_{2}$ are descent statistics: permutation statistics that depend only on the descent composition. We apply this method to a number of descent statistics, including the descent number, the peak number, the left peak number, the number of up-down runs, and the major index. Perhaps surprisingly, in many cases the polynomial giving the joint distribution of $\operatorname{st}_{1}$ and $\operatorname{ist}_{2}$ can be expressed as a simple sum involving products of the polynomials giving the (individual) distributions of $\operatorname{st}_{1}$ and $\operatorname{st}_{2}$. Our work leads to a rederivation of Stanley's generating function for doubly alternating permutations, as well as several conjectures concerning real-rootedness and $\gamma$-positivity.