Understanding the nature of compressible fluctuations in a broad range of turbulent plasmas, from the intracluster medium to the solar wind, has been an active field of research in the past decades. Theoretical frameworks for weakly compressible MHD turbulence in an inhomogeneous background magnetic field predict a linear scaling of the normalized mass density fluctuation ($\delta \rho / \rho_0$), as a function of the turbulent Mach number ($\mathcal{M}_t$), $\delta \rho / \rho_0 \propto \mathcal{M}_t$. However, so far the scaling relation has been tested only using moderate to low plasma beta ($\beta \lesssim 1$) solar wind observational data where the compressibility is weak $\delta \rho / \rho_0 \sim 0.1$. Here, we combine NASA's Magnetospheric Multiscale Mission data in Earth's magnetosheath, where $\beta \sim 10$ is high, and $\beta \sim 1$ highly-compressible magnetohydrodynamic turbulence simulations at unprecedented resolutions. Both show that $\delta \rho / \rho_0 \propto \mathcal{M}_t$ holds across a broad range of $\delta \rho / \rho_0$, $\mathcal{M}_t$ and $\beta$, demonstrating that $\delta \rho / \rho_0 \propto \mathcal{M}_t$ is a robust compressible turbulence relation, going beyond the asymptotics of the weakly compressible theory. We discuss the findings in the context of understanding the nature of strongly compressible turbulent fluctuations and the driving parameter in astrophysical and space plasmas.