A trick for playing a multivariate integral

For suitable functions \(f\) the integral $$ \int_0^{\infty} \int_0^{\infty} \ldots \int_0^{\infty} \frac{dx_1 dx_2 \cdots dx_n}{(f(x_1) + f(x_2) + \ldots + f (x_n))^{k+1}} $$ is equal to $$ \frac{1}{k!} \int_0^{\infty} t^k \left( \int_0^{\infty} e^{- t f(x)} dx \right)^n dt . $$ To see why this is true, start with the last integral and write the \(n\)-th power of the inner integral as a multiple integral over variables \(x_1, \ldots, x_n\). Then change order of integration and evaluate $$ \int_0^{\infty} t^k e^{- t (f(x_1) + \ldots + f(x_n))} dt = \frac{\Gamma(k+1)}{(f(x_1) + \ldots + f(x_n))^{k+1}} $$ which follows from the integral representation $$ \Gamma(z) = \int_0^{\infty} t^{z-1} e^{- t} dt $$ of the gamma function. Note that we don't need to assume that \(k\) is an integer.

I came across this trick while browsing the Journal of Experimental Mathematics and reading the article A Proof of a Recurrence for Bessel Moments. In this article Jonathan M. Borwein and Bruno Salvy are interested in the case where \(f\) is the hyperbolic cosine and the integral $$ K_0(t) = \int_0^{\infty} e^{- t \cosh(x)} dx $$ is a modified Bessel function.