We present a new Bayesian method for reconstructing the parent distributions of X-ray spectral parameters of active galactic nuclei (AGN) in large surveys. The method uses the probability distribution function (PDF) of posteriors obtained by fitting a consistent physical model to each object with a Bayesian method. The PDFs are often broadly distributed and may present systematic biases, such that naive point estimators or even some standard parametric modeling are not sufficient to reconstruct the parent population without obvious bias. Our method uses a transfer function computed from a large realistic simulation with the same selection as in the actual sample to redistribute the stacked PDF and then forward-fit a nonparametric model to it in a Bayesian way, so that the biases in the PDFs are properly taken into account. In this way, we are able to accurately reconstruct the parent distributions. We apply our spectral fitting and population inference methods to the XMM-COSMOS survey as a pilot study. For the 819 AGN detected in the COSMOS field, 663 (81%) of which have spectroscopic redshifts (spec-z) and the others high-quality photometric redshifts (photo-z), we find prominent bi-modality with widely separated peaks in the distribution of the absorbing hydrogen column density (NH) and an indication that absorbed AGN have harder photon indices. A clear decreasing trend of the absorbed AGN fraction versus the intrinsic 2–10 keV luminosity is observed, but there is no clear evolution in the absorbed fraction with redshift. Our method is designed to be readily applicable to large AGN samples such as the XXL survey, and eventually eROSITA.