On the moments of the moments of the characteristic polynomials of random unitary matrices

Denoting by $P_N(A,\theta)=\det(I-Ae^{-i\theta})$ the characteristic polynomial on the unit circle in the complex plane of an $N\times N$ random unitary matrix $A$, we calculate the $k$th moment, defined with respect to an average over $A\in U(N)$, of the random variable corresponding to the $2\beta$th moment of $P_N(A,\theta)$ with respect to the uniform measure $\frac{d\theta}{2\pi}$, for all $k, \beta\in\mathbb{N}$ . These moments of moments have played an important role in recent investigations of the extreme value statistics of characteristic polynomials and their connections with log-correlated Gaussian fields. Our approach is based on a new combinatorial representation of the moments using the theory of symmetric functions, and an analysis of a second representation in terms of multiple contour integrals. Our main result is that the moments of moments are polynomials in $N$ of degree $k^2\beta^2-k+1$. This resolves a conjecture of Fyodorov \& Keating~\cite{fyodorov14} concerning the scaling of the moments with $N$ as $N\rightarrow\infty$, for $k,\beta\in\mathbb{N}$. Indeed, it goes further in that we give a method for computing these polynomials explicitly and obtain a general formula for the leading coefficient.