The study of fluctuations of particle multiplicities in relativistic heavy-ion reactions has drawn much attention in recent years, because they have been proposed as a probe for underlying dynamics and possible formation of quark-gluon plasma. Thus it is of uttermost importance to describe the baseline of statistical fluctuations in the hadron gas phase in a correct way. We performed a comprehensive study of multiplicity distributions in the full ideal hadron-resonance gas in different ensembles, namely grand canonical, canonical, and microcanonical, by using two different methods: Asymptotic expansions and full Monte Carlo simulations. The method based on asymptotic expansion allows a quick numerical calculation of dispersions in the hadron gas with three conserved charges at the primary hadron level, while the Monte Carlo simulation is suitable for studying the effect of resonance decays. Even though mean multiplicities converge to the same values, major differences in fluctuations for these ensembles persist in the thermodynamic limit, as pointed out in recent studies. We observe that this difference is ultimately related to the nonadditivity of the variances in the ensembles with exact conservation of extensive quantities.