  Infinite-dimensional, holomorphic functions have been studied in detail over
the last several decades, due to their relevance to parametric differential
equations and computational uncertainty quantification. The approximation of
such functions from finitely many samples is of particular interest, due to the
practical importance of constructing surrogate models to complex mathematical
models of physical processes. In a previous work, [5] we studied the
approximation of so-called Banach-valued,
$(\boldsymbol{b},\varepsilon)$-holomorphic functions on the
infinite-dimensional hypercube $[-1,1]^{\mathbb{N}}$ from $m$ (potentially
adaptive) samples. In particular, we derived lower bounds for the adaptive
$m$-widths for classes of such functions, which showed that certain algebraic
rates of the form $m^{1/2-1/p}$ are the best possible regardless of the
sampling-recovery pair. In this work, we continue this investigation by
focusing on the practical case where the samples are pointwise evaluations
drawn identically and independently from a probability measure. Specifically,
for Hilbert-valued $(\boldsymbol{b},\varepsilon)$-holomorphic functions, we
show that the same rates can be achieved (up to a small polylogarithmic or
algebraic factor) for essentially arbitrary tensor-product Jacobi
(ultraspherical) measures. Our reconstruction maps are based on least squares
and compressed sensing procedures using the corresponding orthonormal Jacobi
polynomials. In doing so, we strengthen and generalize past work that has
derived weaker nonuniform guarantees for the uniform and Chebyshev measures
(and corresponding polynomials) only. We also extend various best $s$-term
polynomial approximation error bounds to arbitrary Jacobi polynomial
expansions. Overall, we demonstrate that i.i.d.\ pointwise samples are
near-optimal for the recovery of infinite-dimensional, holomorphic functions.
