GAME: GAlaxy Machine learning for Emission lines
Graziano Ucci Andrea Ferrara Andrea Pallottini Simona Gallerani


firstpage–lastpage
On the importance of the nonequilibrium ionization of  and  and 
the line-of-sight in solar surges
    [
    December 30, 2023
==================================================================================================





We present an updated, optimized version of GAME (GAlaxy Machine learning for Emission lines), a code designed to infer key interstellar medium physical properties from emission line intensities of UV/optical/far infrared galaxy spectra. The improvements concern: (a) an enlarged spectral library including Pop III stars; (b) the inclusion of spectral noise in the training procedure, and (c) an accurate evaluation of uncertainties. We extensively validate the optimized code and compare its performance against empirical methods and other available emission line codes (pyqz and HII-CHI-mistry) on a sample of 62 SDSS stacked galaxy spectra and 75 observed HII regions. Very good agreement is found for metallicity. However, ionization parameters derived by GAME tend to be higher. We show that this is due to the use of too limited libraries in the other codes. The main advantages of GAME are the simultaneous use of all the measured spectral lines, and the extremely short computational times. We finally discuss the code potential and limitations.




galaxies: ISM – methods: data analysis – ISM: general, lines and bands, HII regions, PDR






§ INTRODUCTION

Understanding the structure and physical properties of the InterStellar Medium (ISM) in galaxies, especially at high redshift, is one of the major drivers of galaxy formation studies. Measurements of key properties as gas density, column density, metallicity, ionization parameter, Habing flux, rely on galaxy spectra obtained through the most advanced telescopes (both earth-based and spaceborne) and, in particular, on their emission lines <cit.>. However, finding diagnostics that are free of significant systematic uncertainties remains an unsolved problem <cit.>.

Today, extensive spectroscopic studies of the redshifted UV and optical emission lines of distant galaxies (z ≳ 2) are limited to relatively small (with respect to the local Universe) samples of star-forming galaxies <cit.>. The situation is even worse for higher redshifts (z≳ 6) where emission lines measurements are limited to few very bright galaxies <cit.> and gravitationally lensed objects <cit.>. 
One of the main issues is that for such distant galaxies, the number of UV and optical emission lines available is almost always limited to a couple if not to a single extremely bright line. High-quality rest-frame UV and optical spectra including also fainter lines of high-z galaxies, have to await for new facilities, such as the James Webb Space Telescope (JWST) and the Extremely Large Telescope (ELT). These instruments will revolutionize the current research, both in terms of quality and amount of data (estimated production rate is petabyte/year ). 

This latter aspect, in particular, will bring the astrophysical community into an era where Machine Learning (ML) algorithms and Big Data Analytics (BDA) architectures will become fundamental tools in the data mining process. This is already the case for local observations where, e.g., Integral Field Units (IFUs) are already able to provide observations of local galaxies containing tens of thousands of spaxels <cit.>. The development of new methods will be therefore crucial especially in the analysis of galaxy spectra which, combining the efforts from different instruments, will include faint lines arising from very extended wavelengths ranges (i.e. from UV to FIR rest-frame).

The usual emission lines diagnostics will still represent extremely useful tools, albeit they will be likely unable to reach firm conclusions as they rely on a very limited set of lines. These diagnostics are in fact based on UV/optical line ratios <cit.>. This is mostly because UV/optical Recombination Lines (RLs) emitted by the lightest elements (Hydrogen and Helium) are the strongest feature in galaxy spectra. The analogue ones emitted by heavy element can be extremely weak (up to 10^3-4 times fainter than the Hβ 4861 Å line), and thus very difficult to measure in very distant objects. The brightest metal lines are instead the Collisionally Excited Lines (CELs), corresponding to forbidden transitions. 

Some recent works have started to explore new diagnostics based also on FIR lines <cit.>. Facilities such as the Atacama Large Millimeter/submillimeter Array (ALMA) will allow to construct catalogues containing not only information on the UV/optical/near-IR part of the spectrum but also on the far-Infrared (FIR) lines (i.e. [CII] λ157 μm, [OI] λ63 μm, [NII] λ122 μm, [OIII] λ88 μm). Such FIR lines would perfectly complement UV/optical/near-IR lines. However, to fully exploit the information contained in the combined line set, new methods and libraries of synthetic observations are needed.

In this context, several methods in the literature have been developed that widely use ionized gas to infer galaxy physical properties from emission line intensities <cit.>. A common method to measure abundances is based on the determination of the electron temperature. The electron temperature T_e can be calculated from auroral to nebular emission-line ratios such as R_O3 = [OIII] λ4959+5007 / [OIII] λ4363 <cit.> <cit.>. Alternatively, one can use the ratio of RLs of ions which show a weak dependence on T_e and the electron density n_e <cit.>. Given that RLs and auroral lines can be extremely weak in faint, distant or high metallicity objects, other methods which use CELs and Balmer lines have been devised to compute the ISM physical properties. These are referred to as the Strong Emission Line (SEL) methods, and they are based on the comparison of theoretical spectra from a grid of photoionization models <cit.>, or on empirical calibrations obtained for samples for which a previously derived metallicity estimate via the electronic temperature method exists <cit.>. In addition to these, there are codes capable to infer the abundance and ionization parameter of HII regions: IZI <cit.>, pyqz <cit.>, HII-CHI-mistry <cit.>, BOND <cit.>.

Within our group, we have developed GAME (GAlaxy Machine learning for Emission lines), a new fast method to reconstruct the physical properties of the ISM by using all the information represented by the emission lines intensities present in the whole available spectrum <cit.>.  represented a sort of feasibility study: our primary objective was to verify whether ML techniques can be used to predict the physical properties of galaxies. Using the AdaBoost method, we verified that the emission line intensities can effectively provide information on the state of the gas (especially for metallicity). 

As a next step, we present here the current (updated and optimized) version of the GAME code (Section <ref>, for the workflow of the code see also Appendix <ref>) that is based on a new library of photoionization models (50,000 synthetic spectra). In addition to the technical improvements we present also a strategy to deal with uncertainties for a ML-based code. We further implemented an approach to include noise in the library during the ML training phase, in order to apply GAME to real spectroscopic observations. Another key result concerns the treatment of the emission line degeneracy with physical properties. This is discussed in the framework of the application of GAME to study the ISM of star-forming galaxies (Section <ref>). In Section <ref> we also test the performances of GAME comparing it to other methods/codes and against a sample of HII regions with available abundances determinations. We finally discuss the potential and limitations of GAME in Section <ref>.



§ GAME


GAME is a code that, by using as input spectral emission line intensities, infers key galaxy ISM physical properties (see  for full details). It is based on a Supervised Machine Learning algorithm called AdaBoost with Decision Trees as base learner trained with a large library of synthetic spectra. 

To generate each spectrum, we ran the photoionization code CLOUDY v13.03 <cit.>, using as input the quadruplet (n, N_H, U, Z) where n is the total hydrogen density, N_H the column density, U the ionization parameter, and Z the metallicity, We assume an oxygen abundance 12+log( O/H) = 8.69, and solar abundance ratios for all the elements[Such assumption can be relaxed, but this would require to build a dedicated library. We plan to explore the effect of peculiar abundances in future work.] <cit.>. The library covers a large range of physical conditions found in the ISM, as reported in Table <ref>. The total wavelength coverage of the synthetic spectra in the library ranges from the Lyα (1216 Å) wavelength up to 1 mm. However, GAME is able to deal with any subset of emission lines within the input spectra.

A Decision Tree <cit.> recursively partitions the data with respect to the input feature space in branches first (i.e. the branches correspond to different regions of the input feature space), and then into an increasing number of leaves. Because it is possible to improve the power of many base learners into ensemble learning methods <cit.>, we can combine many Decision Trees to make a forest. A common way to produce a forest of Decision Trees is the algorithm called Adaptive Boosting or Adaboost <cit.>. AdaBoost improves the performance of a base learner, by accounting for the elements in the training set which have large prediction errors. We refer the interested reader to <cit.> and <cit.> for a detailed description of the AdaBoost algorithm[Both the implementations of AdaBoost and Decision Trees within GAME are included in the scikit-learn Python package <cit.>, <http://scikit-learn.org/.>].

By running GAME, it is possible to infer four default labels: (n, N_H, U, Z). Besides these default labels (i.e. the physical properties used to generate the CLOUDY photoionization models, see Sec. <ref>), in the GAME library the user can find additional labels, as for example the radius of the cloud r, the visual extinction in magnitudes A_V, or the FUV (6 - 13.6 eV) flux in Habing units G/G_0[G_0 = 1.6× 10^-3 erg  s^-1  cm^-2 <cit.>]. In this paper, we have chosen as additional label G/G_0: the GAME output is therefore a set of values for n, N_H, U, Z, G/G_0.



 §.§ Model improvements


In this Section, we describe some new important improvements introduced here with respect to the original version of the code presented in . These concern (a) the build-up of the spectral library, (b) the inclusion of noise in the training procedure, and (c) the evaluation of the uncertainties; they are described in the following. More technical and detailed materials can be found in the Appendixes.





  §.§.§ Library of synthetic spectra

Concerning the library, the main improvements in the current version are:



  * The library now contains 50,000 synthetic spectra, i.e. it is ∼ 65% larger than that used in ;

  * We added SEDs of Population III stars generated via the YGGDRASIL code <cit.>. We adopted the <cit.> model with a zero-metallicity population and a Kroupa IMF <cit.> in the interval 0.1-100  M_⊙ based on a rescaled single stellar population from <cit.>. For the star formation history we have chosen an instantaneous burst with age set to 2 Myr;

  * In addition to the graphite and silicate dust grains, we have added Polycyclic Aromatic Hydrocarbons (PAH) as these particles considerably affect EUV extinction, and are an important heating source especially in neutral regions. These contributions are modelled following <cit.>. For the effect of stochastic heating we refer to the work of <cit.>. A power-law distribution of PAH is assumed with 10 size bins <cit.>. The dust-to-gas ratio has been linearly scaled with metallicity. 




  §.§.§ Training library with noise


The library discussed so far is made of purely theoretical models. Observed spectra contain noise which must be taken into account. Our approach is to include within the ML training phase a library containing noisy models. To this aim, we have generated a new library (100,000 models) made by two parts: the first is the original library of 50,000 photoionization models. The second is the same library to which gaussian random noise has been added to the lines with an amplitude equal to 10% of the line intensity[Higher noise values would lead to completely noise-dominated models: given the large degeneracy in the emission lines intensities (see Sec. <ref>), GAME would then not be able to reliably discriminate among different models.]. The sum of these two libraries represents therefore our final training dataset. We test GAME against synthetic noisy spectra in Sec. <ref>.



  §.§.§ Uncertainties on the inferred physical properties


The final improvement concerns a robust estimate of the uncertainties associated with the inferred physical properties. The method works as follows. For each input spectrum, we construct N modified versions of it by adding to the line intensities gaussian noise. For each input line intensity I with associated error e, the code extracts N new intensities i from the following Gaussian distribution:


    P(i) = 1/√(2π e^2)exp[-(i - I)^2/2 e^2] .


For upper limits instead, the code generates a new line by taking a random number uniformly distributed between zero and the upper limit.

With this procedure the code generates multiple individual new observations of each spectrum. Then, for each input spectrum we obtain N determinations of the quintuplet of physical properties (n, N_H, G/G_0, U, Z). By default the code gives as final output for each spectrum the average, the median and the standard deviation of these N values. Optionally, the code can return all the N determinations of the physical properties, that can be subsequently combined into a probability distribution function.





§ VALIDATION OF THE CODE




 §.§ Noise in observed spectra

We now analyse how noise on the line intensities can affect the determination of the physical properties.

We start from the emission lines of a synthetic spectrum for which the true physical properties, i.e. those used to generate it, are known. The effect of noise can be mimicked by perturbing each of the synthetic line intensities around its original value. If we then apply GAME to the perturbed intensities, it is possible to assess how the inferred physical properties vary as a function of the noise amplitude. This can be done by using in the ML training phase both the original library (50,000 models) and the one including noise (Sec. <ref>). To perform this analysis we apply the following steps:



  * generate two synthetic spectra: model (a) and model (b) (see Table <ref> for details);

  * choose a set of emission lines to use for the calculation (reported in Table <ref>);

  * choose a set of 20 values of noise percentages (n_i = 1%, 2%, 3%, ..., 20%) and compute 50 different realizations of these two spectra for each n_i value: for each realization we added to the emission line intensities a gaussian random value between 0 and n_i percent of the intensity of the line itself;

  * run GAME to infer the values of the physical properties on both the original library (50,000 photoionization models) and the noisy library (100,000 models, see Section <ref>).






Results are shown in Figure <ref>, where we report the inferred physical properties for models (a) and (b). When the code is trained without noise – first and third columns for models (a) and (b), respectively – adding noise to the spectrum to a level as low as n_i = 4% can lead to differences between the true and inferred properties up to 2 orders of magnitude. Noisier spectra yield even larger differences.

Second and fourth columns of Figure <ref> show that by including the noise in the training procedure reduces these differences in the determination of the physical properties. Interestingly, we also see that metallicity determinations are quite robust. In fact, even adding noise at 20% level the inferred metallicity is within a factor of 2 from the true value, confirming previous conclusions in .

Although the library has been constructed with a noise percentage up to 10%, training the ML algorithm with this library allows GAME to return outputs consistent with the true values even if applied to spectra with greater noise. In fact, as can be seen from Fig. <ref>, the mean of the inferred values agrees quite well with the true values up to noise level as large as 20%.

Using noiseless libraries with low SNR spectra can lead to wrong determinations of the intrinsic physical properties (see Fig. <ref>). This issue is even more severe for weak lines such as [O I] λ6300 or [Ar III] λ7135 (see Table <ref>) which can be fundamental for the determination of the physical properties (i.e. they could have a very high feature importance, see Appendix <ref>).



 §.§ Emission lines degeneracy


Measurements of ISM physical properties are generally based on the comparison between observations and empirically calibrated line ratio or synthetic spectra obtained from photoionization models.







Although empirical calibrations using the electronic temperature are preferable (with respect to theoretical calibrations) because they are based on a quantity directly inferred from observables <cit.>, they also suffer from some limitations <cit.>. The major amongst these is that calibrations often use galaxy samples that do not properly cover all the physical properties space. Hence, empirical calibrations obtained from a sample of low excitation HII regions could give unreliable results when applied to global galaxy spectra <cit.>. 

Comparisons with photoionization models are usually performed changing the metallicity and the ionization parameter but they are limited to a small range of other ISM physical properties (i.e. density) if not only to a single value <cit.>. This is problematic, as the ISM density distribution has a dynamic range that can easily span several dex <cit.>.

For this reason we produced an extended library of physically motivated theoretical models (with 50,000 models described in Section <ref>), whose purpose is to cover the large range of physical conditions found in the ISM. In this Sec. we will show that, with such a library, the emission lines arising from different combinations of the input physical properties (density, column density, ionization parameter and metallicity) are extremely degenerate even if the variation of the physical parameters is at the percent level.

We considered 300 photoionization models with the following physical properties:


    -1.8 <    log(n /  cm^-3) < -1.3
    
    17.3 <    log(N_H /  cm^-2) < 17.8
    
    -1.5 <    log(U) < -1.0
    
    0.3 <    Z / Z_⊙ < 0.6


We then compute the following line intensities: [OIII] λ3727, HeII λ4686, [OIII] λ5007, Hα, [NII] λ6584, [S II] λ6717. In Fig. <ref> we report the ratio of the line intensities for the i-th model and the first one within this set of 300 models. This Fig. shows that, although the physical parameters considered in the models vary within small ranges, still they result into intensity variations of several orders of magnitudes especially for [OIII] λ3727, [NII] λ6584 and [S II] λ6717. 

In Fig. <ref>, in a similar way, we show the variation for some emission line ratios (based on common used emission line diagnostics):


    O_3 N_2   = ([ OIII]λ 5007/ Hβ)/([ NII]λ6584/ Hα) ,
     R_23   = ([ OII]λ 3727 + [ OIII]λ 4959 + [ OIII]λ 5007)/ Hβ ,
     N_2   = [ NII]λ 6584/ Hα ,
     O_32   = [ OIII]λ 5007/[ OII]λ 3727 .


These emission line ratios suffer also from intensity variations up to 2.5 dex. What emerges is that not only the intensity of the lines but also their ratios seem to be affected by a non negligible degeneration. Interestingly, this is not true in the case of R_23, which remains approximately constant, at least in the small range of physical properties considered in this Section. However, it must be noticed that in a larger metallicity range (-2 ≲ log(Z/Z_⊙) ≲ 0.5) R_23 is not monotonically dependent on Z <cit.>. Moreover, one of the fundamental advantages of using GAME, lies in the fact that instead of inferring one single physical property at one time (i.e. metallicity in the case of R_23), it can retrieve simultaneously n, U, N_H, G/G_0, Z for an extended range of physical properties (see Table <ref>).

To fully appreciate this aspect over the entire range of optical wavelengths, using the set of photoionization models listed in Tab. <ref>, we show in Fig. <ref> how a small change in the physical properties is mirrored into large line intensity variations. Panel 1 of Fig. <ref> shows the spectrum for model (a); the remaining panels show the spectral differences with the other models (b)-(f). Although the variation of metallicity between model (a) and (b) is < 0.008 dex, varying simultaneously n, N_H and U by 0.03, 0.08 and 0.07 dex, induces large differences in the resulting spectrum.

GAME is based on a very large library and uses all the information carried by the spectral lines on the physical properties. This approach, allowed by the ML implementation, could overcome the degeneracy better than model fitting techniques, typically based on χ^2 methods, i.e. the minimization of the euclidean distance between models:


    D(m,l) = √(∑_i=1^n (m_i - l_i)^2) ,

where m_i is the value of the emission line intensities of a given model, and l_i the corresponding i-th value contained in the library. In the following we show that minimizing such function does not necessarily lead to a correct determination of the physical properties.





In Fig. <ref> we compare the inferred values of the physical properties using the two approaches, GAME vs. D(m,l) minimization (eq. <ref>), using the set of emission lines in Table <ref>. We perform the test on 10% of the noiseless library (GAME uses the remaining 90% as the training dataset), including 5,000 models. On this reduced set we both infer the physical values with GAME and by distance minimization.

Trying to recover a spectrum that is as similar as possible to the input spectrum (minimizing the distance) leads to worse results with respect to GAME. In Fig. <ref> it is evident that the standard deviation of the logarithm of the ratio between the predicted and true values in the case of GAME is a factor of two smaller than the D(m,l) minimization approach. The similarity between two spectra, given the extreme degeneration (see Sec. <ref>), does not mean necessarily a good correlation with their physical properties. This was expected since other ML techniques as k-Nearest Neighbour (with k=1 in the case of equation <ref>) work well when the input space has a low dimensionality, which does not apply to our problem.



§ COMPARISON WITH OTHER METHODS




 §.§ Overview


The most prominent feature of GAME is that it exploits the full information encoded in a spectrum. Instead of using small, pre-selected subsets of emission line ratios, GAME can use an arbitrary number of lines to infer the ISM physical properties. Additionally, its usage is not limited to a specific ISM phase (i.e. HII regions).

Another key advantage of the ML implementation is that, once trained, it requires a very short computational time (see also Appendix <ref>). This is a crucial feature in view of applications to modern IFU observations with ∼ 100,000 spaxels, each one with a substantial number (> 10) of observed lines per spectrum.



 §.§ Quantitative comparison: galaxies


Here we quantitatively compare GAME with other empirical calibrations and results from different codes. To perform this comparison we use a sample of 62 stacked galaxy spectra. These spectra are obtained from SDSS star forming galaxies in 0.027 < z < 0.25, and are fully described in <cit.>. The spectra are stacked in bins of 0.1 dex according to the log values of their [OII]λ3727/Hβ and [OIII]λ5007/Hβ ratios. For each of the stacked spectra <cit.> derived the corresponding metallicity, using a set of new and self-consistent empirical calibrations. We compare such determination with the results of GAME and the results from two widely used emission line codes, pyqz and HII-CHI-mistry.

pyqz[<http://fpavogt.github.io/pyqz/index.html>] <cit.> is a public Python code which uses abundance- and excitation-sensitive line ratios to define a plane in which the oxygen abundance and ionization parameter can be determined by interpolating a grid of photo-ionization models to match the observed line ratios.

HII-CHI-mistry[<http://www.iaa.es/ epm/HII-CHI-mistry.html>] <cit.> is publicly-available Python code. This code takes the extinction-corrected emission line fluxes and, based on a χ^2 minimization on a photoionization models grid, determines chemical-abundances (O/H, N/O) and ionization parameter. In this work we use the version 3.0 of the code, dealing with input line uncertainties via a Monte Carlo approach.



As evident from the comparison shown in Fig. <ref>, the empirical method tends to slightly overestimate metallicity with respect to all codes. Besides, the mean offset and the standard deviation between GAME predictions and those from other methods are small (less than 0.3 dex). Overall, the agreement on metallicity is good.



We extend the comparison among the three codes to the ionization parameter (Fig. <ref>). Interestingly, the GAME results markedly differ from those obtained with the other two codes. This can be due to the fact that pyqz and HII-CHI-mistry only contain models with a narrower U range (i.e. -3.5 < log U_pyqz < -1.5, and -4.0 < log U_HII-CHI-mistry < -1.5). To verify this point, we ran GAME with a reduced library containing ionization parameter values in the range -4.0 < log(U) < -1.5. Fig. <ref> clearly shows that with the reduced library GAME infers U values that are in agreement with those from the other two codes, e.g. the mean of the inferred values are respectively μ(pyqz) = -2.67, μ(GAME) = -2.76, μ(HII-CHI-mistry) = -2.88. Hence it is important to consider U values that are larger than usually assumed.

Values of U up to 10 are expected in line-emitting galaxies. Let us consider the relation between the ionizing photon flux and the star formation rate <cit.>:


    Q(H) = 2.46 × 10^53(SFR/M_⊙yr^-1)  s^-1


For a compact galaxy with a SFR = 1 M_⊙ yr^-1, n = 1 cm^-3, and radius 1 kpc, we obtain log(U) ∼ 0.6, consistently with the upper limit found with GAME (see Fig. <ref>). As discussed in , the spectrum emerging from a galaxy is weighted by the column densities along the line of sight. If diffuse phases are dominant in the build-up of the final spectrum, their low density pushes U ∝ n^-1 towards large values. 

Note that although the full library contains spectra with -4.0 < log(U) < 3 (see Table <ref>), GAME does not infer values log(U) > 1. This implies that this bound is set by the physical conditions, and it is independent of the library extension.

We finally note that U is the physical parameter most affected by uncertainties. The bootstrap routine now included in the code should however significantly mitigate the problem: although the inferred values for U cover a large range, their PDF obtained with the bootstrap is highly peaked (see Appendix <ref> and Fig. <ref>).



 §.§ Quantitative comparison: HII regions


As an additional comparison, we apply GAME to a sample of 75 observed HII regions with available chemical abundance determinations. We choose HII regions for which a large number (N>13) of high quality (typical errors <10%) emission lines is available. Our final sample is composed by:



  * 66 HII regions located in 21 dwarf irregular galaxies observed by <cit.>. For 25 HII regions, oxygen abundances are obtained via direct detection of emission lines tracing the electron temperature; for the remaining 41 HII regions, abundances are inferred through SEL methods. The emission lines (where available) used as input for GAME are: [OII] λλ3727+3729, [NeIII] λ3869, [OIII] λλ4959+5007, [OI] λ6300, [SIII] λ6312, Hα, [NII] λλ6548+6584, HeI λ6678, [SII] λλ6717+6731, [ArIII] λ7135 (Table 3 of that paper);


  * 9 HII regions located towards the Galactic anti-centre observed by <cit.> for which more than 60 emission lines measurements are available (Table A.1 of that paper).


Results are shown in Fig. <ref>, where we compare the abundance determinations reported in the cited works with the GAME results. The metallicity values inferred by GAME (Z/Z_⊙) are converted into oxygen using the relation 12+log(O/H) = 8.69. The overall agreement between the two methods is good. The scatter does not seem to correlate with metallicity, and the dispersion is of ≃ 0.2 - 0.3 dex. The mean offset (between the results inferred with GAME and the published results reported in Fig. <ref>) and its standard deviation are respectively 0.08 and 0.29. Differences among the results of different methods can arise because of several reasons. First of all, the oxygen abundance we infer from the metallicity assumes solar abundances. Moreover, while GAME adopts input emission lines that are not de-reddened, the reddening corrections adopted by <cit.> and <cit.> (i.e. n_e = 100 cm^-3, Case B recombination, different extinction curves) may have a significant impact on the physical properties determination.



The main advantage of GAME is that it does not use any assumption on the gas physical properties (e.g. density, temperature) defining the emission line spectrum we are looking at. The code is able to recover the physical properties (in this case metallicity) by extracting them from a library containing a vast collection of physical conditions of the ISM. We also note that GAME does not even use the information that the sample refers to HII regions: it simply searches the library, trains itself, and gives the best predicted physical conditions.



§ SUMMARY AND DISCUSSION


We presented an updated and optimized implementation of GAME , a code designed to infer ISM physical properties from emission line spectra. The code is based on a ML algorithm (AdaBoost with Decision Trees as base learner) to calculate density (n), column density (N_H), ionization parameter (U), metallicity (Z) and FUV flux in the Habing band (G/G_0) given an arbitrary set of emission line flux measurements (or upper limits) with their uncertainties.

GAME is extremely reliable particularly for the metallicity determination: the 5-fold cross-validation score with the set of emission lines reported in Table <ref> is higher than 0.95. Although some properties, as metallicity, are easily and robustly recovered also from noisy spectra with few emission lines, other physical properties require higher quality data with more lines measurements. We have shown that for a given set of emission lines if the cross-validation score for the metallicity is 0.95, it might happen that lower scores are obtained for n and N_H (both ∼ 0.8) and U (∼ 0.7).

Another key point in our analysis is that the emission lines intensities are highly degenerate. A small variation of the physical properties leads to large changes in the emission line intensities ratios. The ML approach used here can overcome this issue much better than classical fitting methods based on χ^2 minimization. Noticeably, such important result can be achieved also when the spectra include noise as in real observations.

We have compared GAME with methods based on empirical calibrations <cit.> and other codes (pyqz and HII-CHI-mistry) by considering a sample of 62 stacked spectra from SDSS galaxies <cit.>. While a very good agreement has been obtained in terms of the metallicity determination for the considered sample, we find discrepancies in the derived values of the ionization parameter. We discuss possible reasons for such disagreement in Sec. <ref>.

Finally, we have tested our code on a sample of 75 HII regions with direct method and Strong Emission Lines (SEL) abundances determinations <cit.> to study how GAME can recover these values. We found that the oxygen abundances determinations are in good agreement with those inferred with GAME with a typical scatter around 0.2 - 0.3 dex. The applications of GAME are not only limited to HII regions but the code can deal equally well with different phases of the ISM, including the molecular one. This is because the underlying library covers the largest possible range of physical conditions characterizing the ISM. Furthermore, GAME offers the possibility to use an arbitrary set of emission lines which span a wavelength range from the Lyα one (1216 Å) to 1 mm. These features allow the user to infer the physical properties of phases ranging from the Hot Ionized Medium (HIM) to dense molecular cores.

One of the main limitations of GAME (and other methods/codes) relies on the possible presence of different ISM phases/gradients along the line of sight contributing to the same spectrum. This introduces a complexity that cannot be managed at the present time. We preliminarily discussed this issue in . The main result there was that in such conditions the returned physical parameters are biased towards the phase with largest gas column density. To make progress, we plan to investigate this issue in more details using emission line spectra generated from high-resolution galaxy simulations. Simulated galaxies and their spectra offer the advantage that the physical conditions of the gas shaping the observed spectra is precisely known. This will allow us to (a) devise more stringent reliability tests for GAME, and (b) understand how to maximize the information retrieval from spectra arising from multiphase lines of sight (for this issue we also refer the reader to Section 5.3 of ).

Is is nowadays possible to obtain spatially resolved spectra of galaxies. The information coming from different regions of a galaxy requires the development of new methods in order to obtain the physical conditions of the different phases of the ISM. Large libraries and robust algorithms will be crucial in the analysis of galaxy spectra which include faint lines arising from extended wavelength ranges. Combining the UV/optical/IR/FIR information from the same object will be the next step towards a better understanding of the internal structure of distant galaxies.



§ ACKNOWLEDGEMENTS

We thank M. Curti for providing the galaxy spectra used in our analysis and the anonymous referee for constructive insights. We also thank B. Greig, N. Gillet, C. Behrens, and L. Vallini for useful discussions and comments. AF acknowledges support from the ERC Advanced Grant INTERSTELLAR H2020/740120

mnras







§ WORKFLOW OF THE CODE


Fig. <ref> shows GAME workflow diagram. The diagram presents also an example where 4 input spectra with a maximum of 4 emission lines are used to infer physical properties [We recall that the inferred physical properties are 5, but only 4 are used to generate the Cloudy models: n, N_H, U, Z.]. We underline that GAME can deal with missing lines (third spectrum) or upper limits (fourth). We will refer to this example throughout the Section. The input spectra and the physical property values shown are purely for illustration. 





 §.§ Determination of unique models routine

As the AdaBoost ML algorithm does not deal with missing values, the approach we followed is to simply get rid of emission lines for which a measurement (or at least an upper limit) is not available. The first routine therefore computes the number of unique models (i.e. unique combinations of input emission lines) from the input file. In the example in Fig. <ref> there are two different unique models: (I) one including all the 4 input lines, and (II) another including all but the [OII] line. Each unique model requires a separate  ML training phase. 



 §.§ Accuracy evaluation routine

The goal of any ML model, is to learn from the training set and generalize for unseen data. The accuracy is a measure of the performance on unseen data. To quantify the accuracy, we split the spectral library into a training (90%) and testing (10%) dataset. The accuracy evaluation routine trains the ML algorithm on the training dataset, and predicts values for the testing one. The predicted values are then compared to the actual values used to generate the testing spectra.

In addition, we define a variable k that splits the whole library into k equal parts called folds. The routine then performs the so-called k-fold cross-validation in which the ML is trained on k-1 folds and tested on the remaining fold. 



For each k fold, let p̅_k be the mean of the inferred physical properties:


    p̅_k = 1/n∑_i=1^n p_i,k


where n = (length of library) / k. We then define the two following quantities:


    T_k = ∑_i=1^n (p_i,k - p̅_k)^2 ,



    R_k = ∑_i=1^n (p_i,k - t_i,k)^2 ,

where t_i,k are the true values. For each k and physical property the routine computes the coefficient of determination,


    R^2_k = 1 - R_k/T_k.


The final output is the mean score R^2 =⟨ R^2(k)⟩ and its standard deviation. If the predicted values perfectly fit the true data, R^2=1. Thus R^2 represents a measure of the accuracy of the method.   

Figure <ref> shows the predicted/true ratio distribution for metallicity values using the set of emission lines in Table <ref>. About 90% of the models are within a factor of 2 from the true values, yielding R^2>0.95 (5-fold cross-validation). This is a remarkable performance as it has been obtained with a relatively small number of emission lines. 



 §.§ Machine Learning Routine

This routine trains the AdaBoost ML on the whole library (the one including noise: 100,000 models) and predicts the values for the physical properties. For the example in Fig. <ref> the ML routine is called only twice, one for each of the two unique models (I) and (II). Since there are 5 physical properties to be determined (n, N_H, G/G_0, U, Z), the training procedure is composed by 5 different sub-trainings. We therefore use for the AdaBoost algorithm each label separately and independently.



 §.§ Bootstrap routine

To understand the effects of noise in observed spectra we have devised the following procedure. For each input library spectrum this routine builds N modified versions of it by adding gaussian noise (see Sec. <ref>) to the line intensities. 



We started from a synthetic spectrum featuring the emission lines in Table <ref> (intrinsic values). Then we add gaussian noise with r.m.s. amplitude equal to 15% of the line intensity (error values in Table <ref>).
Fig. <ref> shows the distributions for the inferred physical properties of N = 10,000 realizations[We have verified that the distributions do not change appreciably for larger N values.] for this synthetic input spectrum.



Metallicity is the most accurately-inferred property showing the narrowest distribution (maximum width ∼ 0.2 dex). The inferred and true Z values differ only by 0.02 dex. Although the distributions of the other physical properties are broader, the inferred mean value of, e.g. U, differs only by 0.22 dex from the true one.



§ FEATURE IMPORTANCE

In most cases input features do not contribute equally to predict the outputs. In many situations a certain fraction of the features could be also effectively irrelevant.

The notion of feature importance measures the frequency[Here for frequency we mean the ratio between the number of occurrences of a feature as split point and the total number of split points in the Decision Tree.] of the usage of each feature as split point in a tree, and the resulting purity[In order to decide which feature to split on at each step in building the tree, the algorithm chooses the split that results in the purest daughter nodes. The purity is usually quantified by the reduction of the mean squared error or the entropy.] of the branches determined by a feature separation. Feature importances are defined for the specific ML algorithm used, and the specific set of input lines.

For individual Decision Trees the feature importance interpretation is relatively easy because it is possible to simply look at the tree structure. In fact, a tree intrinsically performs a feature selection by selecting appropriate split points. For the AdaBoost method, which could include several trees, the interpretation becomes more difficult as feature importance is defined as the mean over all trees.

GAME determines the feature importances as real numbers in [0,1] (the feature importances are normalized so that they sum up to 1). The higher the value, the more important is the contribution of the feature to the accuracy. The implementation in the code uses the Gini importance <cit.>. [ More details on the features importance applied to the AdaBoost technique can be found in <cit.>.] Gini importance, sometimes called also mean decrease impurity, is defined as the total decrease in node impurity (weighted by the probability of reaching that node) averaged over all trees of the ensemble <cit.>. 



As an example, we show in Fig. <ref> the feature importance for a subset of emission lines. For all the physical properties considered, the Hα line turns out to be the least important feature. Perhaps surprisingly, we find that the metallicity determination is most affected by the Hβ line. Its importance is ∼ 2 times larger than that of the other emission lines.

From the above, it should be clear that feature importance notion is algorithm-dependent, and does not represent a direct one-to-one connection with a specific physical property. This is because the feature importance of a given line results from the combination of photoionization physics and the architecture of the Decision Trees. For this reason, a purely physical interpretation of such indicator might be misleading.

Nevertheless, the feature importance provides a first indication of the optimal set of lines to be considered as an input for GAME to recover a physical property of interest. For example, lines with high feature importance, if present in the observed spectrum, should be included in the ML analysis as they improve the accuracy of the results. This criterion applies only for lines with a high signal-to-noise ratio, as otherwise the associated noise would limit the benefit of their inclusion. 



§ EXECUTION TIME

The total wall-clock time (on a CPU Intel Core 2 QUAD CPU Q9550/2.83GHz) for running GAME on the synthetic spectrum of Section <ref> with the full library of 100,000 models that includes noise, is around 14 minutes. The sole execution time of the routines that involve the ML algorithm is instead about 13 minutes. These performances have been obtained setting the number of k-folds for the cross-validation score equal to 5, and the total number of trees to 50.

We stress that almost all the running time is used up by the training procedure. If the observed spectra all contain the same lines (i.e. only one unique model) the most time-consuming phase needs to be performed only once. Thus, the total running time will be essentially independent of the number of spectra. 



§ ACCURACY 

We finally discuss how the accuracy of the results depends on the number of lines used and their feature importance. 



We start by showing (Fig. <ref>) the feature importance of the 17 lines listed in Table <ref> for the metallicity label. As in Appendix <ref>, the Hα line turns out to be the least important feature. The 6 most informative features are lines connected to H and He transitions (i.e. Hβ, Hδ, HeI) and 3 metals lines (i.e. [OII], [NeIII], and [ArIII]). Note also that the two most informative lines, Hδ and [ArIII], sum up to >40% of the total feature importance.

Next, we have performed two GAME runs. In run A we analysed 16 different models, with an increasing number of input line intensities (from l=2 to l=17) ordered by their wavelength. Run B is identical, but the lines are added sequentially according to their feature importance (Fig. <ref>).











The results of the metallicity in run A are reported in Figs. <ref> and <ref>. Although it is technically possible to construct Decision Trees even with only two lines, a huge scatter is present in the first two panels of Fig. <ref>. Also noticeable is a thick horizontal line in the distribution, i.e. in many cases GAME predicts the same value. This implies that the code is not effectively learning from data. As a consequence the PDF of the inferred metallicity values for a small number of lines (Fig. <ref>) is essentially flat.

Figs. <ref> and <ref> illustrate[The careful reader might have noticed that the standard deviation, σ, and kurtosis, k, of the PDF with 17 input lines is slightly different from those given in Fig. <ref>. This is because changing the order of the input features in the boosted Decision Trees may result in slightly different Trees.]  the results of run B. Contrary to the previous case, GAME starts to learn already with just two lines, and the scattering remains always smaller than in run A case for any number of input lines. Selecting the input lines according to their feature importance produces clear advantages: using the most important 6 lines produces an accuracy improvement of a factor 2.25 over the wavelength-ordered addition. 

In both runs, the fraction of data outliers drastically drops as l is increased from l=2 to l=10. Beyond that value the standard deviation flattens out to a small value, σ≈ 0.1. To further check this conclusion, We have performed the same study also for all the labels considered here (n, N_H, U, Z, G/G_0). The trend described above is indeed common to all these labels, as it can be appreciated from Fig. <ref>, where we plot the standard deviation as a function of the number of lines used.

In the example of Fig. <ref>, we use 4 input emission lines to infer 5 physical properties, out of which 4 are independent. In the context of traditional statistics, the number of degrees of freedom of the problem is dof = # of input emission lines - # of physical independent properties = 0.

In this Sec. however, we have inferred physical properties using less then 4 emission lines (Figs. <ref> and <ref>). ML and traditional statistics are based on different approaches. ML does not assume any parametric function for describing the data. The method is instead empirical: the outputs or even the applicability of the model are checked a posteriori. Thus, it is not straightforward to extend the dof definition to a ML algorithm. GAME can build trees even with a negative number of dof; it just requires n ≥ 2 lines. A posteriori we can see that only when the code uses more than n ≳ 6 lines the results are statistically reliable (Figs. <ref> and <ref>).





