We investigate the existence of solutions to viscous ergodic Mean Field Games systems in bounded domains with Neumann boundary conditions and local, possibly aggregative couplings. In particular we exploit the associated variational structure and search for constrained minimizers of a suitable functional. Depending on the growth of the coupling, we detect the existence of global minimizers in the mass subcritical and critical case, and of local minimizers in the mass supercritical case, notably up to the Sobolev critical case.