To understand and interpret the observed spectral energy distributions (SEDs) of starbursts, theoretical or semiempirical SED models are necessary. Yet, while they are well founded in theory, independent verification and calibration of these models, including the exploration of possible degeneracies between their parameters, are rarely made. As a consequence, a robust fitting method that leads to unique and reproducible results has been lacking. Here we introduce a novel approach based on Bayesian analysis to fit the Spitzer-Infrared Spectrometer spectra of starbursts using the SED models proposed by Groves et al.. We demonstrate its capabilities and verify the agreement between the derived best-fit parameters and actual physical conditions by modeling the nearby, well-studied, giant H ii region 30 Doradus in the LMC. The derived physical parameters, such as cluster mass, cluster age, interstellar medium pressure, and covering fraction of photodissociation regions, are representative of the 30 Doradus region. The inclusion of the emission lines in the modeling is crucial to break degeneracies. We investigate the limitations and uncertainties by modeling subregions, which are dominated by single components, within 30 Doradus. A remarkable result for 30 Doradus in particular is a considerable contribution to its mid-infrared spectrum from hot (≈300 K) dust. The demonstrated success of our approach will allow us to derive the physical conditions in more distant, spatially unresolved starbursts.