We put established Krylov subspace methods and Gauss quadrature rules to new use by generalizing the class of matrix function-based centrality measures from single-layered to multiplex networks. Our approach relies on the supra-adjacency matrix as the network representation, which has already been used to generalize eigenvector centrality to temporal and multiplex networks. We discuss the cases of unweighted and weighted as well as undirected and directed multiplex networks and present numerical studies on the convergence of the respective methods, which typically requires only few Krylov subspace iterations. The focus of the numerical experiments is put on urban public transport networks.