A random matrix representation is proposed for the two-dimensional (2D) Coulomb gas at inverse temperature β. For 2×2 matrices with Gaussian distribution we analytically compute the nearest-neighbor spacing distribution of complex eigenvalues in radial distance. Because it does not provide such a good approximation as the Wigner surmise in 1D, we introduce an effective β_{eff}(β) in our analytic formula that describes the spacing obtained numerically from the 2D Coulomb gas well for small values of β. It reproduces the 2D Poisson distribution at β=0 exactly, that is valid for a large particle number. The surmise is used to fit data in two examples, from open quantum spin chains and ecology. The spacing distributions of complex symmetric and complex quaternion self-dual ensembles of non-Hermitian random matrices, that are only known numerically, are very well fitted by noninteger values β=1.4 and β=2.6 from a 2D Coulomb gas, respectively. These two ensembles have been suggested as the only two symmetry classes, where the 2D bulk statistics is different from the Ginibre ensemble.