The distribution of masses for neutron stars is analyzed using the Bayesian statistical inference, evaluating the likelihood of proposed gaussian peaks by using fifty-four measured points obtained in a variety of systems. The results strongly suggest
the existence of a bimodal distribution of the masses, with the first peak around $1.37 {M_{odot}}$, and a much wider second peak at $1.73 {M_{odot}}$. The results support earlier views related to the different evolutionary histories of the members for the first two peaks, which produces a natural separation (even if no attempt to label the systems has been made here), and argues against the single-mass scale viewpoint. The bimodal distribution can also accommodate the recent findings of $sim M_{odot}$ masses quite naturally. Finally, we explore the existence of a subgroup around $1.25 {M_{odot}}$, finding weak, if any, evidence for it. This recently claimed low-mass subgroup, possibly related to $O-Mg-Ne$ core collapse events, has a monotonically decreasing likelihood and does not stand out clearly from the rest of the sample.
