We present, to our knowledge, the first Bayesian extraction of baryon number susceptibilities of QCD matter at finite baryon density from heavy-ion collision data on proton number cumulants. The framework embeds arbitrary equation-of-state susceptibilities $χ_n^B$ into realistic hydrodynamic particlization hypersurfaces through maximum-entropy freeze-out, maps the resulting (anti)baryon fluctuations onto protons, applies the experimental kinematic acceptance, and accounts for exact baryon number conservation. Applying the framework to measurements of (net-)proton number fluctuations in 0–5% central Au-Au collisions from the RHIC Beam Energy Scan in the collider mode, we extract, at each collision energy, the second-order susceptibility normalized by the hadron resonance gas value, $χ_{2}^B/\barχ_{2}^B$, and the higher-order susceptibility ratios $χ_{3}^B/χ_{1}^B$ and $χ_{4}^B/χ_{2}^B$. We obtain tight constraints on $χ_{2}^B$, with extracted values in quantitative agreement with lattice QCD based estimates along the chemical freeze-out line for $μ_B \lesssim 300$~MeV. At larger $μ_B$, the extracted values indicate an enhancement of baryon number fluctuations relative to the noncritical lattice-based extrapolation and HRG baseline considered here. The third- and fourth-order susceptibilities are only weakly constrained. In particular, we find that the observed nonmonotonic collision-energy dependence of the proton factorial cumulant ratio $\hat{C}{3}/\hat{C}{1}$ can be described without irreducible three- or four-baryon correlations. This behavior emerges from the interplay between the energy dependence of $χ_{2}^B$ and exact baryon number conservation.