We study the relaxation of a highly collisional, ultracold but nondegenerate gas of polar molecules. Confined within a harmonic trap, the gas is subject to fluid-gaseous coupled dynamics that lead to a breakdown of first-order hydrodynamics. An attempt to treat these higher-order hydrodynamic effects was previously made with a Gaussian ansatz and coarse-graining model parameter [R. R. W. Wang & J. L. Bohn, Phys. Rev. A 108, 013322 (2023)], leading to an approximate set of equations for a few collective observables accessible to experiments. Here we present substantially improved reduced-order models for these same observables, admissible beyond previous parameter regimes, discovered directly from particle simulations using the WSINDy algorithm (Weak-form Sparse Identification of Nonlinear Dynamics). The interpretable nature of the learning algorithm enables estimation of previously unknown physical quantities and discovery of model terms with candidate physical mechanisms, revealing new physics in mixed collisional regimes. Our approach constitutes a general framework for data-driven model identification leveraging known physics.