  We present a new method for inferring the metallicity (Z) and ionization
parameter (q) of HII regions and star-forming galaxies using strong nebular
emission lines (SEL). We use Bayesian inference to derive the joint and
marginalized posterior probability density functions for Z and q given a set of
observed line fluxes and an input photo-ionization model. Our approach allows
the use of arbitrary sets of SELs and the inclusion of flux upper limits. The
method provides a self-consistent way of determining the physical conditions of
ionized nebulae that is not tied to the arbitrary choice of a particular SEL
diagnostic and uses all the available information. Unlike theoretically
calibrated SEL diagnostics the method is flexible and not tied to a particular
photo-ionization model. We describe our algorithm, validate it against other
methods, and present a tool that implements it called IZI. Using a sample of
nearby extra-galactic HII regions we assess the performance of commonly used
SEL abundance diagnostics. We also use a sample of 22 local HII regions having
both direct and recombination line (RL) oxygen abundance measurements in the
literature to study discrepancies in the abundance scale between different
methods. We find that oxygen abundances derived through Bayesian inference
using currently available photo-ionization models in the literature can be in
good (~30%) agreement with RL abundances, although some models perform
significantly better than others. We also confirm that abundances measured
using the direct method are typically 0.2 dex lower than both RL and
photo-ionization model based abundances.
