Extending the wavenumber-explicit analysis of [Chen & Qiu, J. Comput. Appl. Math. 309 (2017)], we analyze the $L^2$-convergence of a least squares method for the Helmholtz equation with wavenumber $k$. For domains with an analytic boundary, we obtain improved rates in the mesh size $h$ and the polynomial degree $p$ under the scale resolution condition that $hk/p$ is sufficiently small and $p / \log k$ is sufficiently large.