We prepare and study strongly interacting two-dimensional Bose gases in the superfluid, the classical Berezinskii-Kosterlitz-Thouless (BKT) transition, and the vacuum-to-superfluid quantum critical regimes. A wide range of the two-body interaction strength 0.05 < g < 3 is covered by tuning the scattering length and by loading the sample into an optical lattice. Based on the equations of state measurements, we extract the coupling constants as well as critical thermodynamic quantities in different regimes. In the superfluid and the BKT transition regimes, the extracted coupling constants show significant down-shifts from the mean-field and perturbation calculations when g approaches or exceeds one. In the BKT and the quantum critical regimes, all measured thermodynamic quantities show logarithmic dependence on the interaction strength, a tendency confirmed by the extended classical-field and renormalization calculations.PACS numbers: 51.30.+i, 67.25.dj, 64.70.Tg, 37.10.Jk Two-dimensional (2D) Bose gases are an intriguing system to study the interplay between quantum statistics, fluctuations, and interaction. For noninteracting bosons in 2D, fluctuations prevail at finite temperatures and Bose-Einstein condensation occurs only at zero temperature. The presence of interaction can drastically change the picture. With repulsive interactions, fluctuations are reduced and superfluidity emerges at finite temperature via the Berezenskii-Kosterliz-Thouless (BKT) mechanism [1,2]. Interacting Bose gases in two dimensions and BKT physics have been actively investigated in many condensed matter experiments [3][4][5][6][7]. In cold atoms, the BKT transition and the suppression of fluctuations are observed based on 2D gases in the weak interaction regimes [8][9][10][11].Extensive theoretical research on 2D Bose systems addresses the role of interactions in the superfluid phase [12][13][14][15][16][17] and near the BKT critical point [18,19]. In the weak interaction regime, the classical φ 4 field theory [18,19] predicts the logarithmic corrections to the critical chemical potential µ c = k B T (g/π) ln(13.2/g) and the critical density n c = λ −2 dB ln(380/g) for small two-body interaction strength g < 0.2. Here k B T is the thermal energy and λ dB is the thermal de Broglie wavelength. The classical-field predictions are consistent with weakly interacting 2D gas experiments [9][10][11]20]. The upper blue shaded area is the superfluid regime, and the red boundary corresponds to the BKT transition regime. The black dashed lineμ = 0 indicates where we evaluate the density and pressure for a vacuum-to-superfuid quantum critical gas. The inset compares the equations of state of a 2D gas and a 2D lattice gas with an almost identical g ≈ 0.4. µ MF = 2 gn/m logarithmically [13]. Here, m is the mass of the boson, n is the density, and 2π is the Planck constant. Defining the superfluid coupling constant as G = m/( 2 κ), where κ = ∂n/∂µ is the compressibility, we can summarize the perturbation expansion result of G as [12]