All classes of Active Galactic Nuclei (AGN) are fundamentally powered by accretion of gas onto a supermassive black hole. The process converts the potential energy of the infalling matter to X-ray and ultraviolet (UV) radiation, releasing up to several 1012 solar luminosities. Observations show that the accreting "central engines" in AGN are surrounded by dusty matter. The dust occupies a "torus" around the AGN which is comprised of discrete clumps. If the AGN radiation is propagating through the torus on its way to an observer, it will be heavily re-processed by the dust, i.e. converted from UV to infrared (IR) wavelengths. Much of the information about the input radiation is lost in this conversion process while an imprint of the dusty torus is left in the released IR photons. Our group was the first to formulate a consistent treatment of radiative transfer in a clumpy medium an important improvement over simpler models with smooth dust distributions previously used by researchers. Our code CLUMPY computes spectral energy distributions (SED) for any set of model parameters values. Fitting these models to observed AGN SEDs allows us to determine important quantities, such as the torus size, the spatial distribution of clumps, the torus covering factor, or the intrinsic AGN luminosity. Detailed modeling also permits us to study the complex behavior of certain spectral features. IR radiative transfer introduces degeneracies to the solution space: different parameter values can yield similar SEDs. The geometry of the torus further exacerbates the problem. Knowing the amount of parameter degeneracy present in our models is important for quantifying the confidence in data fits. When matching the models to observed SEDs we must employ modern statistical methods. In my research I use Bayesian statistics to determine the likely ranges of parameter values. I have developed all tools required for fitting observed SEDs with our large model database: the latest implementation of CLUMPY, the fit algorithms, the Markov Chain Monte Carlo sampler, and the Bayesian estimator. In collaboration with observing groups we have applied our methods to a multitude of real-life AGN.