  We prove that the grand canonical Gibbs state of an interacting quantum Bose
gas converges to the Gibbs measure of a nonlinear Schr\"odinger equation in the
mean-field limit, where the density of the gas becomes large and the
interaction strength is proportional to the inverse density. Our results hold
in dimensions $d \leq 3$. For $d > 1$ the Gibbs measure is supported on
distributions of negative regularity and we have to renormalize the
interaction. More precisely, we prove the convergence of the relative partition
function and of the reduced density matrices in the $L^r$-norm with optimal
exponent $r$. Moreover, we prove the convergence in the $L^\infty$-norm of
Wick-ordered reduced density matrices, which allows us to control correlations
of Wick-ordered particle densities as well as the asymptotic distribution of
the particle number. Our proof is based on a functional integral representation
of the grand canonical Gibbs state, in which convergence to the mean-field
limit follows formally from an infinite-dimensional stationary phase argument
for ill-defined non-Gaussian measures. We make this argument rigorous by
introducing a white-noise-type auxiliary field, through which the functional
integral is expressed in terms of propagators of heat equations driven by
time-dependent periodic random potentials and can, in turn, be expressed as a
gas of interacting Brownian loops and paths. When the gas is confined by an
external trapping potential, we control the decay of the reduced density
matrices using excursion probabilities of Brownian bridges.
