The nonparametric probabilistic method (NPM) is a physics-based machine learning approach for model-form (MF) uncertainty quantification (UQ), model updating, and digital twinning. It extracts from reference data information not captured by a real-time computational model and infuses it into a “hyperparameterized,” stochastic version of this model by solving an inverse statistical problem formulated using a loss function and a few hyperparameters. The loss function is designed such that the mean value and statistical fluctuations of some quantities of interest predicted using the real-time stochastic model match target values obtained from the reference data. NPM’s performance hinges upon the efficient minimization of the loss function, which is unfortunately stochastic and nonconvex and thus prone to getting the optimization procedure trapped in suboptimal local minima. The latter scenario is exacerbated when the reference data is scarce. The paper addresses these issues by adopting the squared quadratic Wasserstein distance as the measure of dissimilarity between two different sets of data and by reformulating NPM’s inverse statistical problem as a multimodal data-assimilation problem. The potential of the resulting enhanced NPM for MF-UQ, model updating, and digital twinning is demonstrated using numerical simulations relevant to applications in structural dynamics, including structural health monitoring.