Issue 
A&A
Volume 657, January 2022



Article Number  A19  
Number of page(s)  20  
Section  Extragalactic astronomy  
DOI  https://doi.org/10.1051/00046361/202140414  
Published online  21 December 2021 
Scaling relations and baryonic cycling in local starforming galaxies
III. Outflows, effective yields, and metal loading factors
^{1}
INAF – Osservatorio Astronomico di Capodimonte, Salita Moiariello 16, 80131 Napoli, Italy
email: crescenzo.tortora@inaf.it
^{2}
INAF – Osservatorio Astrofisico di Arcetri, Largo Enrico Fermi 5, 50125 Firenze, Italy
^{3}
European Southern Observatory, KarlSchwarzschildStraße 2, 85748 Garching, Germany
^{4}
Observatoire de Genève, Université de Genève, 51 Ch. des Maillettes, 1290 Versoix, Switzerland
Received:
25
January
2021
Accepted:
4
October
2021
Gas accretion and stellar feedback processes link metal content, star formation, and gas and stellar mass (and the potential depth) in starforming galaxies. Constraining this hypersurface has been challenging because of the need for measurements of H I and H_{2} gas masses spanning a broad parameter space. A recent step forward has been achieved through the Metallicity And Gas for Mass Assembly sample of local starforming galaxies, which consists of homogeneously determined parameters and a significant quantity of dwarf galaxies, with stellar masses as low as ∼10^{5} − 10^{6} M_{⊙}. Here, in the third paper of a series, we adopt a standard galactic chemical evolution model, with which we can quantify stellardriven outflows. In particular, we constrain the difference between the massloading in accretion and outflows and the wind metalloading factor. The resulting model reproduces very well the local mass–metallicity relation, and the observed trends of metallicity with gas fraction. Although the difference in mass loading between accreted and expelled gas is extremely difficult to constrain, we find indications that, on average, the amount of gas acquired through accretion is roughly the same as the gas lost through bulk stellar outflows, a condition roughly corresponding to a “gas equilibrium” scenario. In agreement with previous work, the wind metalloading factor shows a steep increase toward lower mass and circular velocity, indicating that lowmass galaxies are more efficient at expelling metals, thus shaping the mass–metallicity relation. Effective yields are found to increase with mass up to an inflection mass threshold, with a mild decline at larger masses and circular velocities. A comparison of our results for metal loading in outflows with the expectations for their mass loading favors momentumdriven winds at low masses, rather than energydriven ones.
Key words: galaxies: star formation / galaxies: ISM / galaxies: fundamental parameters / galaxies: statistics / galaxies: dwarf / evolution
© ESO 2021
1. Introduction
Starforming galaxies evolve and assemble their mass by accreting gas from the intergalactic medium (IGM), converting it in stars, then expelling gas and metals, thus changing the gas metallicity in the interstellar medium (ISM). This baryon cycle is mainly regulated by internal processes^{1} (e.g., Wu et al. 2017), which are responsible for metal and dust production, and the subsequent expulsion of material. Outflows driven by stellar winds and supernovae (SNe) explosions dominate in the lowermass regime, while at higher masses, outflows tend to be powered by feedback from active galactic nuclei (AGN, e.g., Tremonti et al. 2004; Dekel & Birnboim 2006; Dalcanton 2007; Erb 2008; Finlator & Davé 2008; Zahid et al. 2014; De Rossi et al. 2017; LaraLópez et al. 2019).
The efficiency of these processes is predicted to change with galaxy properties, and in particular with the virial mass (M_{vir}) and stellar mass (M_{star}) of the galaxies. These variations lead to the most important relations among the main galaxy parameters driving galaxy evolution, namely the “main sequence” correlation between M_{star} and the starformation rate (SFR, e.g., Brinchmann et al. 2004, the SFMS); the mass–metallicity relation (e.g., Tremonti et al. 2004, the MZR); and the consequent inverse correlation between metallicity and specific SFR (sSFR ≡ SFR/M_{star}, e.g., Mannucci et al. 2010; LaraLópez et al. 2010; Hunt et al. 2012, 2016). Baryonic processes are thought to shape the SFMS and the MZR, in addition to many other correlations among other stellar properties and dark matter: the virialtostellar mass (e.g., Moster et al. 2010) and the baryonic Tully–Fisher (e.g., Lelli et al. 2016a) relations; the correlation between color and stellar population gradients with mass (e.g., Tortora et al. 2010); and the relation of dark matter fraction and mass density slopes with stellar mass (e.g., Tortora et al. 2019).
These scaling relations are not necessarily straight power laws, but can also show curvatures, bending, breaks, and Ushapes, suggesting the existence of typical mass scales which regulate the efficiency of the physical processes involved in baryonic cycling (e.g., Dekel & Birnboim 2006). Kannappan et al. (2013) identified two transitions in galaxy morphology, gas fractions, and fueling regimes. The first transition was considered a “gasrichness” threshold; galaxies below a transition mass of log(M_{star}/M_{⊙}) ∼ 9.5 were found to be “accretiondominated”, characterized by overwhelming gas accretion. The second M_{star} transition corresponds, instead, to a bimodality threshold, roughly indicating a “quenched” regime at log(M_{star}/M_{⊙}) > 10.5, populated by spheroiddominated, gaspoor galaxies.
Both M_{star} thresholds, and the mass regimes they delineate, are really tracing gas content, and its intimate link with starformation activity (Kannappan et al. 2013). To further explore this link, in the first paper of this series (Ginolfi et al. 2020, hereafter Paper I), we constructed a sample of galaxies in the Local Universe with homogeneously derived global measurements of the main parameters that quantify galaxy evolution: namely M_{star}, SFR, metallicity O/H, and gas content including molecular gas mass M_{H2} and atomic gas mass M_{HI}. The sample comprises ∼400 galaxies, and was dubbed Metallicity And Gas for Mass Assembly (MAGMA).
In a second paper (Hunt et al. 2020, hereafter Paper II), we used MAGMA to explore how molecular and atomic gas properties vary as a function of M_{star} and SFR. The main focus of Paper II was on feedback processes, and how they govern scaling relations across the mass spectrum. Essentially three processes play a role in regulating galaxy evolution^{2}: (P1) preventive feedback, the (lack of) availability of cold baryons from the host halo; (P2) inefficiency of the starformation process, the conversion of the available gas into stars; and (P3) ejective feedback, the production of energy and momentum, and the expulsion of material including gas, dust, and metals. In Paper II, we investigated two of the three main feedback mechanisms behind the inefficiency of lowmass halos to form stars (e.g., Behroozi et al. 2013; Moster et al. 2013; Graziani et al. 2015, 2017): namely the availability of cold gas to form stars (P1), and the inefficiency of the starformation process itself (P2).
Here, in the third paper of the series, we examine the remaining mechanism, (P3), the efficiency of metal enrichment and ejective feedback. Understanding the processes that shape galaxy scaling relations with metallicity has proved to be a challenging task. It has been historically important to compare observations with the “closedbox” model, assuming that the total baryon content in the galaxy is fixed and the gas and stellar content are only determined by star formation and the stellar yields (e.g., Pagel & Patchett 1975). However, it is widely recognized that galaxies do not evolve as closed boxes (e.g., Tremonti et al. 2004), and more complex models have considered inflows and outflows of material from the galactic ecosystem (e.g., Larson 1972; Pagel & Patchett 1975; Tinsley 1980; Pagel & Edmunds 1981; Edmunds & Pagel 1984; Edmunds 1990; Tremonti et al. 2004; Dalcanton 2007; Erb 2008; Finlator & Davé 2008; Recchi et al. 2008; Peeples & Shankar 2011; Dayal et al. 2013; Zahid et al. 2014; Peng & Maiolino 2014). More recently, hydrodynamical simulations have tested these predictions (e.g., Ma et al. 2016; De Rossi et al. 2017; Torrey et al. 2019).
Our approach here is mainly empirical, requiring a sample with measurements of M_{star}, SFR, O/H, and gas content, homogeneously determined over a wide range of stellar mass and gas fractions; MAGMA is such a sample. In this paper, Sect. 2 describes briefly the MAGMA sample and how it was assembled. Then, Sect. 3 presents the theoretical formalism for metallicity evolution and the scaling of metallicity with M_{star}, SFR, and gas content. The computation of outflow parameters and metalloading factors for MAGMA galaxies is described in Sect. 4, and their impact on shaping the mass–metallicity relation in Sect. 5. We recast our results in the context of the effective yields in Sect. 6, and discuss the physical implications for wind metallicity and mass loading of the winds in Sect. 7. Our conclusions are drawn in Sect. 8.
2. MAGMA sample
Paper I presented the MAGMA sample, comprising 392 local galaxies having homogeneously derived global quantities including M_{star}, SFR, metallicities [12 + log(O/H)], and gas masses (both atomic, M_{HI}, and molecular, M_{H2}, and total M_{g} = M_{HI} + M_{H2}). Molecular gas masses were obtained by measurements of the CO luminosity, . We assembled this sample by combining a variety of previous surveys at z ∼ 0 with new observations of CO in lowmass galaxies. Only galaxies with robust detections of physical parameters were retained, that is to say we required clear detections of CO and H I. Galaxies were eliminated if they showed optical signatures of an AGN in the parent samples (obtained from lineratio diagnostics, similar to Baldwin et al. 1981). We also eliminated potentially H Ideficient galaxies by applying a cut based on H Ideficiency measurements, where available (e.g., Boselli et al. 2009, 2014a). The MAGMA sample spans more than 5 orders of magnitude in M_{star}, SFR, and M_{g}, and a factor of ∼50 in metallicity [Z, 12 + log(O/H)].
Because of potential systematics that could perturb our results, in Paper I we homogenized the sample by recalculating both M_{star} and SFR in a uniform way. To calculate SFR, we adopted the hybrid formulations of Leroy et al. (2019), based on combinations of GALEX and WISE luminosities. New M_{star} measurements were computed using 3.4 μm (or IRAC 3.6 μm when not available) luminosities from fluxes in the ALLWISE Source Catalogue (Wright et al. 2010) and a luminositydependent M/L from Hunt et al. (2019); the continuum contribution from freefree emission was subtracted from these luminosities before the M_{star} computation. Both SFR and M_{star} assume a Chabrier (2003) IMF.
For metallicity, 12 + log(O/H), we adopted either “direct” calibrations, based on electron temperatures (e.g., Izotov et al. 2007), or the [N II]based strongline calibration by Pettini & Pagel (2004, PP04N2). When PP04N2derived metallicities were not available, they were converted from the original calibration to PP04N2 according to Kewley & Ellison (2008). PP04N2 is the metallicity calibration closest to the direct one (e.g., Andrews & Martini 2013; Hunt et al. 2016), which makes it appropriate for a sample like MAGMA that spans a wide range of M_{star} and SFR.
The molecular gas mass was determined from based mainly on ^{12}CO(1 − 0); only in a few cases did we have to convert from ^{12}CO(2 − 1) to ^{12}CO(1 − 0) using “standard” prescriptions (e.g., Leroy et al. 2009; Schruba et al. 2012). The conversion factor α_{CO} from CO luminosity to molecular gas mass was calculated according to the piecewise relation presented in Paper II, with α_{CO} ∝ Z^{1.55}, for 12 + log(O/H) < 8.69 and (K km s^{−1} pc^{2})^{−1} (see also Saintonge et al. 2011). The impact of a different metallicity dependence is discussed in Paper II. Here we do not include a factor of 1.36 to account for helium. For more details about the parent surveys, as well as parameter determinations and calibrations, quality checks, or the derivation of scaling relations, we refer the reader to Paper I and Paper II.
3. Metallicity evolution and scaling relations: formalism
Much effort has been spent over the last fifty years to understand the theoretical framework behind the baryon exchange cycle and metallicity evolution in galaxies (e.g., Larson 1972; Pagel & Patchett 1975; Tinsley 1980; Pagel & Edmunds 1981; Edmunds & Pagel 1984; Edmunds 1990; Tremonti et al. 2004; Dalcanton 2007; Erb 2008; Finlator & Davé 2008; Recchi et al. 2008; Spitoni et al. 2010; Peeples & Shankar 2011; Dayal et al. 2013; Zahid et al. 2014; Peng & Maiolino 2014). It has been fairly well established that bulk outflows and galactic winds driven by star formation are necessary to explain the MZR and other scaling relations involving metallicity. Here we follow these pioneering works, and explore the behavior of the MZR, SFR, and gas content in the MAGMA galaxies. MAGMA provides the important constraints of gas measurements together with SFR and metallicity, thus enabling, for the first time, a quantification of the mass and metallicity loading in stellar outflows, and their role in shaping the MZR in the Local Universe.
3.1. Galactic chemical evolution
In what follows, we rely heavily on the formulations for galactic chemical evolution by Pagel (2009). We first make the assumption of instantaneous recycling, namely that the chemical enrichment provided by stellar evolution, nucleosynthesis, and expulsion of metals into the ISM take place “instantaneously” relative to the timescales of overall galaxy evolution.
There are two main quantities that this assumption makes timeinvariant: the first is the return fraction R, the mass fraction of a single stellar generation that is returned to the ISM; the second is the (true) stellar yield, y, defined as the mass of an element newly produced and expelled by a generation of stars relative to the mass that remains locked up in longlived stars and compact remnants. The metal mass produced relative to the initial mass of a single stellar generation is then q = α y where α = 1 − R, (or “lockup” fraction), with the implicit assumption that stellar lifetimes can be neglected so that R is a net mass return fraction (see e.g., Vincenzo et al. 2016).
These parameters are governed by the IMF, and in particular the mass fraction of a stellar generation above a certain limit, since these are the main contributors to metal production. We further assume that metals are homogeneously mixed at all times, and that initally there are no stars and no metals, only the initial gas M_{g}(t = 0) with mass M_{i}.
There are four main physical quantities involved in and constrained by chemical evolution: the mass in stars, M_{star}; the starformation rate SFR, here referred to as ψ; the mass in gas, M_{g}; and the mass in metals, M_{Z}, or the metal abundance in the gas, Z_{g} = M_{Z}/M_{g}. The mass of stars that have been born up to time t can be written as:
Consequently, the mass M_{star} remaining in stars and longlived remnants is:
The total baryonic mass in a galactic system is the sum of gas and stars:
where M_{w} is the mass expelled in bulk winds and outflows, and M_{a} is the mass of newly accreted gas. This is essentially a mass continuity equation that defines the mass conservation between infall, outflows, and star formation. Thus, the time derivative of the gas mass can be written as:
The time variation of the mass of metals M_{Z} in the ISM gas depends on ψ, R, the stellar yield q, Ṁ_{a}, and Ṁ_{w} as follows:
where Z_{a} and Z_{w} are the metallicities of the accreted and outflow gas, respectively. The first two terms on the righthand side of Eq. (5) correspond to the addition of metals of ejecta from stars; the third term to the loss to the ISM from astration; and the remaining terms to the metals gained from accretion and lost in bulk stellar outflows from winds and SNe.
Equation (5) can be simplified by two common assumptions: we take Ṁ_{a} and Ṁ_{w} to be directly proportional to ψ (e.g., Erb 2008; Peeples & Shankar 2011; Dayal et al. 2013; Lilly et al. 2013):
where we have introduced the massloading factors η_{a} and η_{w}. It is fairly well established that the latter relation is true, namely that the mass in outflows is proportional to SFR, at least to its integral over time (e.g., Muratov et al. 2015; Christensen et al. 2016). In contrast, the former proportionality is more subject to speculation, although evidence is mounting that gas accretion powers star formation (e.g., Fraternali & Tomassetti 2012; Sánchez Almeida et al. 2014; Tacchella et al. 2016). In any case, this assumption provides the noteable advantage of analytical simplicity; we plan to investigate more completely its ramifications in a future paper. Following Peeples & Shankar (2011), we also introduce the metalloading factors ζ_{a} and ζ_{w} relating Z_{g}, Z_{a}, η_{a}, Z_{w}, and η_{w}:
and
ζ_{w} and ζ_{a} essentially describe the efficiency with which gas accretion and starformationdriven outflows change metal content within the galaxy. Thus, we can rewrite Eqs. (4) and (5) as:
and
Our observations do not directly constrain M_{Z}, but rather Z_{g}, so we need to reformulate Eq. (10):
in order to solve for Ż_{g}:
where we have used Eq. (9) for Ṁ_{g} and Eq. (10) for Ṁ_{Z}.
The necessary ingredients are now in place to derive the solution for the dependence of Z_{g} on the remaining parameters. First, we know that SFR is related to gas mass by:
where ϵ_{⋆} is the inverse depletion time or timescale for star formation. Equations (9) and (13) are coupled, so must be solved simultaneously. Assuming that ϵ_{⋆} does not vary with time, we can thus integrate Eq. (9), to obtain the following solution for the gas mass:
Substituting this solution into Eq. (12) and integrating, we obtain:
where it can be seen that the time dependence is encapsulated in the evolution of the gas (M_{g}/M_{i}, e.g., Recchi et al. 2008). Equation (15) is essentially telling us that galaxies evolve along a hypersurface constrained by M_{star}, M_{g}, gas metallicity, and accretion and wind mass and metal loadings. The galaxy’s position on this hypersurface can change with time as gas is accreted, stars are formed, and mass and metals are ejected in winds. However, at each time, the relation constraining these fundamental variables is approximately invariant.
The ratio of the system gas mass M_{g} and the initial gas mass M_{i} can be written in terms of the baryonic gas mass fraction, μ_{g} = M_{g}/(M_{g} + M_{⋆}):
where we have assumed that the total baryonic mass is described by Eq. (3). This is again an expression of mass conservation, meaning that no gas has been expelled altogether from the galaxy, but rather has been recycled into stars as described in Eq. (13).
Similar solutions for Z_{g} were first formulated by Edmunds & Pagel (1984), and later by Recchi et al. (2008), Erb (2008), Spitoni et al. (2010), Peeples & Shankar (2011), Dayal et al. (2013), and many others. These earlier works mostly assumed that Z_{w}/Z_{g} = 1 in a socalled homogeneous wind (e.g., Pagel 2009), leading to ζ_{w} = η_{w}. This assumption simplifies Eq. (15) placing the dependence of η_{w} only in the exponent of M_{g}/M_{i}. However, as has been shown previously (e.g., Peeples & Shankar 2011) and as we show below for MAGMA, this is almost certainly an oversimplified assumption because the metallicity of the outflowing material Z_{w} is not generally the same as the ISM gas, Z_{g}, especially at low masses. We further discuss this simplification and its implications later on in the paper (see also Appendix C).
3.2. Limitations of the model
First, Eq. (15) is highly sensitive to the value of η_{a} − η_{w} − α in the denominator of the exponent of M_{g}/M_{i}. For convenience in what follows, we define:
When Δη = α, the exponent for M_{g}/M_{i} in Eq. (15) becomes infinite. This equality, Δη ≈ α, is exactly what would be expected for the “equilibrium scenario” for gas and metallicity evolution, in which gas accretion is balanced by star formation and winds (e.g., Davé et al. 2012; Mitra et al. 2015), such that Ṁ_{g} ≈ 0. From Eq. (9), this implies that either ψ ≈ 0 which is not physically relevant here, or Δη ≈ α. In other words, for “equilibrium”, the difference between the accretion η_{a} and wind η_{w} coefficients is balanced by the lockup fraction α. However, this formalism gives for the equilibrium scenario a solution for Z_{g} that is apparently illconditioned^{3}. When Δη = α, Eq. (14) (and Eq. (16)) are irrelevant, as M_{g} = M_{i}.
More generally, the trend of gas growth in Eq. (14) depends on the sign of the exponent Δη − α. This exponent is a critical component of Eq. (15) because of the sensitivity to the sign of the exponent of M_{g}/M_{i}. Where Ṁ_{g} > 0, with Δη − α slightly positive, a finite Z_{g} solution is obtained only when this exponent is negative, and conversely, if Δη − α is slightly negative, then this exponent should be positive.
Equation (15) is also sensitive to the gas fraction μ_{g} through Eq. (16). As also noticed by Recchi et al. (2008), this approach leads to constraints on the possible values of gas fraction μ_{g}. Depending on the value of Δη − α, very small μ_{g} leads to infinite solutions for Z_{g}.
In any case, despite the complexity of possible solutions, Eq. (15) is a useful tool to assess the interplay of gas and metallicity in galaxies. In the next section, we explore this by letting the data guide the possible solutions.
4. Mass and metal loading in MAGMA
MAGMA combines homogeneous measurements of M_{star}, SFR, Z_{g}, and M_{g}; so, for the first time, we can infer from observations the mass loading factors, η_{a} and η_{w}, as well as the metalloading factor ζ_{w}. To do this, we first make another common assumption, namely that the accreted gas is pristine (Z_{a} = 0, ζ_{a} = 0, e.g., Erb 2008; Peeples & Shankar 2011; Dayal et al. 2013; Lilly et al. 2013; Creasey et al. 2015). Although this assumption may not be exactly true for lowredshift galaxy populations, there is considerable observational evidence that the gas accreted by local galaxies is significantly metal poor (e.g., Tosi 1982; Sánchez Almeida et al. 2014; Fraternali 2017).
To apply the formalism described above to the observations, we convert the MAGMA oxygen abundance 12 + log(O/H) to Z_{g} according to:
where the underlying assumption is that oxygen measures the enrichment of a primordial mix of gas comprising approximately 75% hydrogen by mass (e.g., Peeples & Shankar 2011; Chisholm et al. 2018). We adopt the yields and return fractions for the Chabrier (2003) IMF from Vincenzo et al. (2016) using the Romano et al. (2010) stellar yields: (y, R) = (0.037, 0.455). We have purposely not considered a conversion of O/H to total metal abundance Z_{g}, because of the inherent uncertainties in this conversion (e.g., Vincenzo et al. 2016); thus our value for y is the oxygen yield y_{O} from Vincenzo et al. (2016).
From Eq. (15), and fixing ζ_{a} to 0, because of our assumption that inflowing gas is pristine, it can be seen that there are two variables that shape the Z_{g} − M_{star}–SFR–M_{g} hypersurface: Δη and ζ_{w}. With this formulation, it is impossible to determine separately the accretion and outflow massloading factors, η_{a} and η_{w}. Physically, this means that the important driving factor is the difference Δη, that is to say the dominance (or not) of gas accretion over bulk outflows. Below we describe the results from a Bayesian approach for inferring mass and metallicityloading factors in MAGMA.
4.1. Bayesian approach
For each galaxy with four observables (M_{star}, SFR, Z_{g}, M_{g}), we need to determine two parameters, Δη and ζ_{w}. However, only two observables enter into Eq. (15): Z_{g} and μ_{g}. Thus, the problem is severely underdetermined. To overcome this limitation, we have divided MAGMA into M_{star} bins, and assumed that Δη and ζ_{w} are the same for each mass bin. The M_{star} bins were defined to guarantee a sufficient number of galaxies in each bin, but also to adequately emphasize the lowestM_{star} MAGMA galaxies.
To find the bestfit combination of (Δη, ζ_{w}) for each M_{star} bin, we first calculate χ^{2} value for each of N galaxies in the bin^{4}:
where for the jth galaxy is given by Eq. (18) and by Eqs. (15) and (16) for the ith parameter pair (Δη, ζ_{w}); we have assumed a constant error σ = 0.1, similar to the deviation of the MZR fit discussed in Paper I. χ^{2} is computed over the grid of Δη and ζ_{w} parameters, in order to obtain N values (indexed by j in Eq. (19)) of χ^{2} for each parameter pair (Δη, ζ_{w}).
Since only one value of χ^{2} can be associated with each unique parameter pair, a decision must be made for the assignment of this single value for the N galaxies in the mass bin. This can sometimes be accomplished through stacking (e.g., Belfiore et al. 2019), but we cannot do that here, because of the way is defined (Eqs. (15) and (16)). We experimented with several methods; the approach that best captures the sensitivity of the χ^{2} distribution to the variation in the parameters (Δη, ζ_{w}) is to assign the single χ^{2} value to the mean of the lowest quartile in the χ^{2} distribution of the N galaxies in the mass bin. This is a rather arbitrary choice but the numbers of galaxies in each M_{star} bin are sufficient to sample in a meaningful way the galaxies with the best fits. Thus for each mass bin, assigning this single χ^{2} value to each parameter pair, we obtain the distribution of χ^{2} values for that bin.
We then apply the standard Bayesian formulation:
to derive the full posterior probability distribution P(θD) of the parameter vector θ = (Δη, ζ_{w}), given the data vector D = (Z_{g}, μ_{g}). This posterior is proportional to the product of the prior P(θ) on all model parameters (the probability of a given model being obtained without knowledge of the data), and the likelihood P(Dθ) that the data are compatible with a model generated by a particular set of parameters. The data are assumed to be characterized by Gaussian uncertainties, so the likelihood of a given model is proportional to exp(−χ^{2}/2). We assume uniform priors: Δη ∈ [ − 3.4, 3.4] in linearly sampled steps, and log_{10}(ζ_{w}) ∈ [ − 0.1, 2.5] in logarithmic steps of 0.05 dex.
The bestfit parameter vector (Δη, ζ_{w}) is evaluated by constructing the probability density functions (PDFs), weighting each model with the likelihood exp(χ^{2}/2), and normalizing to ensure a total probability of unity. Specifically, this was achieved by marginalization over other parameters:
where Y is either Δη or ζ_{w}, excluding the parameter of interest (ζ_{w} or Δη, respectively).
The results for the different mass bins are shown in Fig. 1. The χ^{2} surfaces show the relation between Δη and ζ_{w}. For each M_{star} bin, the PDF median is shown as a vertical dashed line, and the 1σ confidence intervals (CIs) by a shaded violet region. We have taken the maximumlikelihood estimate (MLE) to be the median of the PDF for ζ_{w} and the mode for Δη. This choice is discussed more fully below.
Fig. 1. Corner plots of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for MAGMA binned by M_{star}. The violet contours correspond to the minimum χ^{2} value. The top and right panels of each corner plot report the probability density distributions for the marginalized parameters; confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF median) is shown by a vertical dashed line. 
Figure 1 shows that Δη is degenerate; low values of χ^{2} can be achieved for a broad range of Δη. However, strikingly, the spread of allowed values for Δη at low M_{star} (in the two upper leftmost panels in Fig. 1) is much broader than for higher stellar masses. In the higher M_{star} bins, the allowed distribution of Δη is truncated to increasingly small values of Δη. This implies that η_{a}, with equal probability, can be significantly larger than η_{w} at low masses, as also shown by the PDFs in each panel. Another clear result that emerges from Fig. 1 is that higher ζ_{w} occur at lower M_{star}.
Figure 1 also shows that while the PDF of ζ_{w} shows a form similar to a Gaussian distribution, and is thus fairly well constrained, the PDF for Δη is not. It peaks at a given value, but then falls off toward lower values gradually, leading to the conclusion that the bestfit value may depend on the interval over which the uniform prior is sampled. To investigate this, we have performed the Bayesian calculations using two different intervals in Δη:Δη ∈ [ − 10, 3.4] and Δη ∈ [0, 3.4], again in linearly sampled steps. These cornerplots are shown in Appendix A, and the bestfit values for Δη and ζ_{w} for Δη ∈ [0, 3.4] are given in the last two columns of Table 1.
Δη and ζ_{w} obtained in bins of mass for Bayesian method.
Results are shown graphically in Fig. 2, where Δη is plotted against V_{circ} in the upper panel and ζ_{w} vs. V_{circ} in the lower one. Because the outflow strengths probably depend on the depth of the potential well, we assess them in terms of the virial velocities, V_{circ}, as proposed by Peeples & Shankar (2011). We have calculated virial velocities in several ways, but ultimately adopt V_{out}, velocities from the baryonic Tully–Fisher relation (BTFR) predicted by APOSTLE/EAGLE simulations from Sales et al. (2017), focused on the lowmass regime. Other alternatives are also viable, such as the powerlaw BTFR from Lelli et al. (2016b), or velocities obtained from the M_{vir} − M_{star} relation via halo abundance matching (Moster et al. 2010). However, at low M_{star}, and low baryonic masses, the inferred virial velocities change significantly according to whether the relation is a straight power law as in Lelli et al. (2016b) or inflected at low masses as in Sales et al. (2017). The simple powerlaw relation gives V_{circ} values that are apparently too low relative to the rest of the sample, while the V_{circ} values by Sales et al. (2017) are more consistent. Moreover, because the lowmass end of our sample is somewhat extreme, there is also some doubt about the applicability of the usual haloabundance matching techniques. In any case, our conclusions are independent of the formulation used to translate the mass scale to the kinematic one.
Fig. 2. Difference of massloading factors Δη (top panel) and metalloading factors ζ_{w} (bottom panel) as a function of virial velocities obtained from the BTFR from Sales et al. (2017). In both panels, filled symbols show the Bayesian results, and in the top panel, open ones illustrate the upper bound of the PDF 1σ confidence intervals (CI). Different symbol types correspond to the interval of Δη in the Bayesian priors, as described in the text; the symbols for Δη ∈ [0, 3.4] and Δη ∈ [ − 10, 3.4] are (arbitrarily) offset along the abscissa for better visibility. The gray regions give the ±1σ uncertainties in the bestfit values. Also shown is the “equilibrium” asymptote, namely Δη = α (see Sect. 3.1). In the top panel, the two estimates of Δη from different priors are exactly coincident. In the bottom panel, the formulations by Peeples & Shankar (2011) for two O/H calibrations, redone here for our value of y, are given by (purple) shortdashed and (blue) longdashed curves for the PP04N2 and D02 calibrations, respectively. Also shown are the metalloading factors computed by Chisholm et al. (2018) from COS data. The text gives more details. 
As shown in Fig. 2, the different intervals used for the Bayesian priors give fairly consistent results both for Δη and for ζ_{w}. The bestfit Bayesian value for Δη is ≈α for all the prior intervals sampled, and there is a clear trend for ζ_{w} to increase toward lower V_{circ} (or stellar mass, as we show later). Over the mass range sampled by MAGMA, ζ_{w} is almost a factor of ten higher at lower M_{star} (V_{circ}) relative to the most massive galaxies.
Also shown in the bottom panel of Fig. 2 are the estimates of ζ_{w} obtained by Chisholm et al. (2018) with Cosmic Origins Spectrograph (COS) observations of ultraviolet absorption lines. By analyzing the underlying stellar continuum, and using the observed line profiles to constrain optical depth and covering fraction, Chisholm et al. (2018) were able to quantify the mass outflow rate, Ṁ_{w}. Then, through photoionization modeling, they obtain column densities of each element that are then used to constrain the metallicity of the outflow Z_{w}. As shown in Chisholm et al. (2018), their results are well fit by the analytical predictions of Peeples & Shankar (2011), with yields in the range 0.007−0.038, which agree with our reference value of y = 0.037.
4.2. Comparison with a closedform approximation
Figure 2 also illustrates the predictions of Peeples & Shankar (2011) using their formulation for the MZR (taken from Kewley & Ellison 2008, and converted to a Chabrier 2003 IMF). Peeples & Shankar (2011) inferred ζ_{w} from the behavior of M_{g} and Z_{g} with M_{star}, using known derivatives derived empirically; we discuss this approach more fully below. Figure 2 shows their formulation relating ζ_{w} and V_{circ}:
where, for consistency, we have recomputed their fits using our higher adopted yield (for y_{O}) of 0.037 taken from Vincenzo et al. (2016) for a Chabrier (2003) IMF. As shown in Fig. 2, the PP04N2 results are in good agreement with MAGMA, while the D02 calibration (Denicoló et al. 2002) is slightly discrepant, and more consistent with Chisholm et al. (2018), despite the higher yield. This illustrates how the computed ζ_{w} values depend on the adopted yields, and underscores the necessity of a selfconsistent framework. We have performed our entire metallicity analysis with lower (and higher) yields than the ones used here, and the results do not change.
Peeples & Shankar (2011) assessed Δη and ζ_{w} based on known derivatives calculated from observed scaling relations of M_{g} and Z_{g} with M_{star}. In some sense, this is complementary to our Bayesian approach, but may also skew results as we show below. Following Pagel (2009), the time dependence in Eq. (9) can be eliminated by noting that Ṁ_{⋆} = α ψ. Thus, by changing variables we obtain:
The same can be done for Eq. (10), and eliminating the time dependence as above gives:
Our observations do not directly constrain M_{Z}, but rather Z_{g}, so we need to formulate dZ_{g}/dM_{⋆} as follows:
The lefthand side is given by Eq. (24), and dM_{g}/dM_{⋆} by Eq. (23), so we can solve for dZ_{g}/dM_{⋆}:
As they should be, Eqs. (23) and (26) are related to Eqs. (9) and (12) by the simple derivative, Ṁ_{⋆} = α ψ. Equations (23) and (26) are slightly different from the ones derived by Peeples & Shankar (2011), because of the different assumptions they made for the stellar enrichment given by Eq. (5).
For a given M_{g} vs. M_{star} relation, we know dM_{g}/dM_{star}, the left side of Eq. (23); for a given MZR, we know dZ_{g}/dM_{star}, the left side of Eq. (26). Thus, in the spirit of Peeples & Shankar (2011), Δη = η_{a} − η_{w} and ζ_{w} can be directly constrained observationally. Inverting Eqs. (23) and (26) gives Δη and ζ_{w} explicitly:
In Eq. (28), we have defined F_{g} ≡ M_{g}/M_{star}. Peeples & Shankar (2011) were the first to derive such a formulation for ζ_{w}; as mentioned in Sect. 3, the formalism here differs slightly from theirs because of the way they dealt with metal production in stellar ejecta.
The problem with Eq. (27) is that for positive d log M_{g}/d log M_{star}, Δη is by definition > 0, because F_{g} is ≥0. In Paper II, the robust leastsquares power–law relation^{5} between M_{g} and M_{star} is found to be:
thus with d log M_{g}/d log M_{star} ≈ 0.7, roughly consistent with earlier works (e.g., Leroy et al. 2008). Because of this positive derivative, to infer Δη by relying on Eq. (27) and d log M_{g}/d log M_{star} would give results that are inherently skewed toward higher values of Δη. Qualitatively, we would expect that as a galaxy evolves, the stellar mass continues to grow, while the gas mass diminishes; such a scenario would give d log M_{g}/d log M_{star} < 0, and possibly Δη < α. Consequently, from a physical point of view, it may be inaccurate or unrealistic to invoke the sample’s overall d log M_{g}/d log M_{star} to calculate Δη. What we require is an estimate of Δη for individual galaxies which may or may not be related to the statistical properties of the parent sample at z ≈ 0.
In any case, Fig. 2 illustrates that ζ_{w} is not particularly sensitive to Δη. The predictions of ζ_{w} from Peeples & Shankar (2011) using d log M_{g}/d log M_{star} are fairly consistent with those we find with the Bayesian approach, even though the former requires Δη > 0 and the latter suggests that Δη ≈ 0.
5. Shaping the mass–metallicity relation with outflows
In the previous section, we estimated Δη and ζ_{w} through a Bayesian approach. Here we parameterize Δη and ζ_{w} as a function of both M_{star} and μ_{g} to predict quantitatively how mass loading Δη and metal loading ζ_{w} shape the MZR.
5.1. Stellar mass, metallicity, and outflows
First, we incorporate M_{star} as the defining variable. Figure 3 (left panel) shows the trends of Bayesianinferred Δη and ζ_{w} with observed Z_{g} as a function of M_{star} (in logarithmic space). The bottom left panel shows the predictions for Z_{g} that are inferred from the functional forms we have found for Δη and ζ_{w} (see Table 2). The middle left panel shows ζ_{w} together with the robust cubic bestfit polynomials with coefficients as given in Table 2. The cubic polynomial has no physical significance; rather our objective is to parameterize Δη and ζ_{w} relative to M_{star} in order to compute Z_{g}.
Fig. 3. Difference of massloading factors Δη and metalloading factors ζ_{w} inferred from the Bayesian analysis, and observed Z_{g} in the bottom panel, as a function of M_{star} on the left, and as a function of μ_{g} on the right. In the two upper panels, the bestfit Bayesian values are shown as filled symbols (circles or squares according to the Δη prior), with error bars corresponding to the ±1σ excursion, which is also illustrated by lightgray shaded regions. The rightmost point (magenta circle) in the top and middle right panels corresponds to the “extra” highμ_{g} bin shown in Fig. 6. In the top panels, a horizontal dashed line illustrates the “equilibrium solution” with Δη = α. In the top left panel, the (purple) dashed curve is a linear powerlaw fit of Δη to log(M_{star}), while the dashed curve in the top right panel shows Eq. (27), with d log M_{g}/d log M_{star} = 0.72, as found in Paper II. The fits to the curves for ζ_{w} in the middle panels are cubic: log(ζ_{w}) as a function of log(M_{star}) in the left panel, and ζ_{w} as a function of F_{g} in the right. In the lower panels, observed MAGMA Z_{g} are shown as small filled circles, color coded by gas fraction, while the Z_{g} medians and 1σ deviations are shown by large (black) filled circles and error bars. The dashed curves in the lower panels correspond to the calculation of Z_{g} based on the fits in the upper panels (as a function of μ_{g} on the left, and log(M_{star}) on the right): the longdashed curves (dark orange) correspond to the Δη ≥ 0 priors, while the shortdashed ones (purple) to the “symmetric” priors. The green curves in the left panels are described in the text, obtained by fitting individual galaxy Δη and ζ_{w} calculated using Eqs. (27) and (28). 
Robust polynomial fits to Bayesian estimates of Δη and ζ_{w}.
The slight differences (∼0.15 dex for M_{star} ≳ 10^{9} M_{⊙}, longdashed orange vs. shortdashed purple curves) in ζ_{w} according to the Δη prior are evident in the middle left panel in Fig. 3 and the lower one. When the constraint on the prior Δη > 0, the Bayesian inference of ζ_{w} gives lower values, and is a better approximation of observed Z_{g} (orange longdashed curve). On the other hand, when the Δη priors include negative values ζ_{w} is larger (see also Fig. 2), but Z_{g} is underestimated (purple shortdashed curve). Since the curves shown in the lower panel are not fits to Z_{g}, but rather the parameterization of Z_{g} as a function of Δη and ζ_{w}, this may be telling us that positive values over negative Δη are preferred, at least in the current evolutionary states of the MAGMA galaxies.
The green curves in the top two left panels of Fig. 3 illustrate the polynomial fits (quadratic for Δη and cubic for log ζ_{w}) that would be obtained had we calculated Δη and ζ_{w} from Eqs. (27) and (28), respectively, using d log Z_{g}/d log M_{⋆} and d log M_{g}/d log M_{⋆} derived from the MAGMA scaling relations. The corresponding green curve in the lower left panel is the approximation of Z_{g} that would be obtained from these fits. It is clear that for log M_{star}/M_{⊙} ≳ 8.5, Z_{g} is equally well estimated by the priors of Δη ≥ 0 and the fits to individual parameters computed with Eqs. (27) and (28). Wind metal loading ζ_{w} estimated with the symmetric prior of Δη falls short (∼0.15 dex) of observed Z_{g} in this mass range. In this mass range, Fig. 3 also shows that Δη is unimportant; the main factor that shapes Z_{g} is ζ_{w}.
On the other hand, the fit to Δη in the top panel (green curve, Eq. (27)) is significantly rising toward lower M_{star}, unlike the Bayesian values for Δη that, with both sets of priors, is roughly constant ≈α. The influence of Δη at low M_{star} (log M_{star}/M_{⊙} ≲ 8.5) is clearly shown by the green curve in the lower left panel of Fig. 3. The falloff of Z_{g} toward low M_{star} is less severe than would be predicted by a flat Δη, and is more consistent with the observations.
5.2. Gas fraction, metallicity, and outflows
In Paper II (see also Eq. (29)), we found that gas content F_{g} grows with decreasing stellar mass (see also e.g., Haynes et al. 2011; Huang et al. 2012; Boselli et al. 2014b; Saintonge et al. 2017; Catinella et al. 2018). Figure 4 shows this comparison explicitly for baryonic gas fraction μ_{g} and M_{star} in MAGMA galaxies. There is a clear trend for higher μ_{g} as M_{star} decreases, but in the very lowest M_{star} regime, below log(M_{star}/M_{⊙}) ≲ 8.5, the scatter is large. In MAGMA, most galaxies in this regime have μ_{g} ≳ 0.5, but there are two exceptions, Sextans A and WLM, at very low μ_{g} despite low global stellar mass. This apparent inconsistency probably arises from underestimates of the total gas content corresponding to the localized regions where CO was detected (e.g., Shi et al. 2015, 2016; Elmegreen et al. 2013; Rubio et al. 2015), and indicates the inherent difficulties of probing the gas, stellar, and metal contents in these extreme dwarf systems. Although the median trend for lessmassive galaxies to be more gas rich is evident down to M_{star} ∼ 3 × 10^{8} M_{⊙}, the large overall scatter makes it virtually impossible to infer trends of Δη and ζ_{w} with μ_{g} from their trends with M_{star}.
Fig. 4. Baryonic gas fraction μ_{g} plotted against (log) M_{star} for MAGMA galaxies. Galaxies are coded by gas fraction, F_{g}, as shown in the vertical color bar. MAGMA medians are shown as a heavy gray line, with ±1σ excursions as the gray region. As noted in the text, the two lowmass galaxies with low μ_{g} are Sextans A and WLM, with probable substantial underestimates of the total gas content. 
Thus, to assess the dependence of Δη and ζ_{w} on μ_{g}, we cannot use the analysis of Sect. 4.1. Consequently, we repeated the Bayesian estimation by binning in μ_{g}, rather than in M_{star}. As before, bin boundaries were chosen such that a similar number of galaxies is found in each of the nine bins; the priors for Δη and ζ_{w} are as before (see Sect. 4.1). The results for the symmetric prior on Δη ∈ [ − 3.4, 3.4] are shown in Fig. 5, where the layout of the figure is the same as in Fig. 1. Here there is a trend for larger Δη excursions at high μ_{g} (see lower panels), similar to what we found for low M_{star}. There is also a tendency for the distribution of ζ_{w} to be broader at high μ_{g}, that is also reflected in the broader range of Δη.
Fig. 5. Corner plots of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for MAGMA binned by μ_{g}. The violet contours correspond to the minimum χ^{2} value. The top and right panels of each corner plot report the probability density distributions for the marginalized parameters; confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF mode for Δη, median for ζ_{w}) is shown by a vertical dashed line. 
At high μ_{g}, there is some evidence for a bimodal distribution of ζ_{w}, as shown in Fig. 6, where we have raised the boundary of the highest μ_{g} bin. In this μ_{g} bin, the rightmost one plotted in Fig. 5, the bestfit (mode) Δη is higher than any. Nevertheless, as can be seen in the lowerright panel, higher values of ζ_{w} are possible for virtually every Δη sampled, so that the dependence of ζ_{w} on Δη is apparently degenerate, as we concluded in Sect. 4.1.
Fig. 6. Corner plot of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for the highest μ_{g} bin MAGMA, here with a higher boundary than in Fig. 5: μ_{g} > 0.75. As in previous figures, the violet contours correspond to the minimum χ^{2} value, and the top and right panels report the probability density distributions for the marginalized parameters. Confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF mode for Δη, median for ζ_{w}) is shown by a vertical dashed line. 
To parameterize Δη and ζ_{w}, here we adopt F_{g} as the independent variable because it is statistically better behaved than μ_{g}, as a result of the simpler dependence on M_{star}. However, it is more convenient to compare quantities to μ_{g} as shown in the right panel of Fig. 3 where Bayesian estimates for Δη, ζ_{w}, and observed Z_{g} are plotted as a function of μ_{g}. The curve in the top right panel is not a fit, but rather the prediction of Eq. (27) using d log M_{g}/d log M_{star} = 0.72 from Eq. (29). The middle right panel shows ζ_{w} together with the robust cubic bestfit polynomial with coefficients as given in Table 2. Metallicity vs. μ_{g} is given in the bottom right panel, where the longdashed curve corresponds not to a fit, but rather to the prediction of Eq. (15) with the functional forms for ζ_{w} (cubic polynomial) and Δη. The difference in ζ_{w} for the two different priors for Δη in the Bayesian analysis are analogous to those for M_{star} (left panel Fig. 3). At low μ_{g}, Z_{g} is better approximated by ζ_{w} inferred from Δη > 0 relative to ζ_{w} from the symmetric prior.
The influences of ζ_{w} and Δη are clearly seen in the trend for Z_{g} with μ_{g}. The steep increase of ζ_{w} and Δη toward high μ_{g} is reflected in the dramatic reduction of metallicity toward high μ_{g}. The shape of the MAGMA Z_{g} curve relative to μ_{g} resembles that found by Erb (2008, see her Fig. 3), for z ∼ 2 galaxies, but with some important differences. The first is that Erb (2008) assume that Z_{w}/Z_{g} = 1, that is to say that ζ_{w} = η_{w}. Appendix C presents the trivial derivation of this result. Under this assumption, Z_{g} in Eq. (15) is essentially governed by η_{a}, with η_{w} entering only into the M_{g}/M_{i} exponent (see also Appendix C). Consequently, there is less variation at low μ_{g}, with the prediction by Erb (2008) assuming a constant asymptotic value for μ_{g} as high as ∼0.35. By considering the change in Z_{w} relative to Z_{g}, as we do here, Z_{g} never reaches a truly asymptotic value, although there is flattening toward lower μ_{g}.
The second difference is the assumption for R (thus α) in the metallicity treatment of Erb (2008). She uses the same metallicity calibration as adopted here (PP04N2), but assumes that α, the fraction of gas retained in stars, is unity, thus neglecting the gas that is expelled in SNe and stellar winds. This would imply a difference −0.34 dex in her and our absolute metallicity scale. An additional factor is the different yield adopted here (y = 0.037) relative to y = 0.019 used by Erb (2008). Considering these two numerical differences, the absolute scale in Erb (2008) would be reduced by ∼0.6 dex. Such a difference is consistent with the galaxies shown by Erb (2008), considering that her sample is at z ∼ 2, while MAGMA is at z ∼ 0. Thus, as observed, lower Z_{g} (∼ − 2.7 dex) at low gas fractions would be expected for her sample relative to ours (e.g., Erb et al. 2006; Maiolino et al. 2008).
5.3. Considerations
Interpretation of our results is not altogether straightforward. On the one hand, the Bayesian analysis shows that Δη as a function of M_{star} tends toward α which would indicate a condition of “gas equilibrium” (e.g., Davé et al. 2012; Mitra et al. 2015). On the other, Δη as a function of μ_{g} does not show this behavior, but rather resembles what would be predicted by Eq. (27). This distinction stems from the highly imperfect correlation of μ_{g} with M_{star} as shown in Fig. 4. Figures 1 and 5 also illustrate that a wide range of Δη is possible (Δη < 0, Δη < α, Δη ≥ α), but that this only weakly affects ζ_{w}.
Our analysis also shows that within a 1σ range, Δη can be < α, over all mass ranges, but possibly more importantly, Δη can be significantly > α toward low M_{star} and high μ_{g} (see Fig. 3). We take this to imply that at low M_{star} and high μ_{g} there is more possibility for increased overall gas accretion (η_{a}) relative to wind outflows (η_{w}). This could be related to the stochastic nature of starformation histories in lowmass galaxies (e.g., Weisz et al. 2012; Madau et al. 2014; Krumholz et al. 2015; Emami et al. 2019), and to the availability of gas.
Figure 3 shows that Δη ≈ α almost independently of M_{star}, thus implying gas equilibrium at all masses; however, Δη increases as a function of μ_{g} roughly as predicted by Eq. (27). Because we would expect that younger, less evolved, galaxies contain more gas, thus higher μ_{g}, this would seem to suggest that Δη could depend on the evolutionary phase of a galaxy. Since we cannot know the starformation histories (SFHs) of the individual galaxies, this suggestion is difficult to test. However, we can fix a simple SFH, such as a falling or rising exponential, and experiment with what such a “toy model” would predict for d log M_{g}/d log M_{star} as a function of time. If we consider that the timescale of star formation τ_{SF} is tightly linked to the gas depletion time τ_{g} (), we can explicitly calculate the growth of M_{star} and M_{g} using Eqs. (2) and (14), respectively.
These toymodel calculations are shown in Appendix B. As expected, d log M_{g}/d log M_{star} depends clearly on Δη; both positive and negative values of d log M_{g}/d log M_{star} are possible depending on the SFH, the age of the galaxy t_{gal}, the timescale for star formation τ_{SF}, and most importantly on Δη. The value of Δη, relative to α in the exponent of Eq. (14), determines whether gas mass is increasing or diminishing with time. Since stellar mass is monotonically increasing, the value of Δη sets the value of d log M_{g}/d log M_{star}. Nevertheless, our results show that Δη is almost impossible to constrain, demonstrating behavior that varies according to whether the independent variable is M_{star} or gas content μ_{g}.
As discussed in Sect. 3.2, when Δη = α, as for the “gasequilibrium” scenario, Eq. (15) is numerically intractable because of the zero denominator in the exponent for M_{g}/M_{i}. This same condition gives Ṁ_{g} = 0 and M_{g} = M_{i} (see Eq. (16)). Because of the relatively long totalgas depletion times, ∼a few Gyr, even without a true equilibrium solution, essentially M_{g}/M_{i} remains unity throughout most of the lifetime of the galaxy (see Eq. (14)).
An important aspect of our formalism is that, in essence, through Eq. (13), we assume that the SFH corresponds also to the time evolution of gas mass. Moreover, because we assume that total system mass is conserved and that the time derivative of accreted gas mass Ṁ_{a} is proportional to SFR through η_{a}, it is possible that the formalism itself requires a rough equality between Δη and α, the lockup fraction.
Despite the degeneracy of Δη and the inability of our formalism to determine unambiguously its behavior, ζ_{w} is fairly well determined (see Figs. 2 and 3). After discussing effective yields in Sect. 6, in the remainder of the paper, we focus on ζ_{w}, and what can be learned from our results about metal and massloading in the stellardriven outflows of starforming galaxies.
6. Effective yields
Above we have quantified a scenario for the relation between Z_{g} and gas content in galaxies that requires significant gas accretion and stellardriven outflows. In this section, we recast our results in the context of the effective yield, as measured against the predictions from a different, “closedbox”, paradigm. A galaxy evolving as a closed box obeys a simple relationship between the metallicity and the gas mass fraction (Searle & Sargent 1972; Pagel & Patchett 1975):
where μ_{g} is the baryonic gas mass fraction and y is defined as the “true” yield (Sect. 3.1). As gas is converted into stars, metal content grows and the gas mass diminishes. If the galaxy evolves as a closed box, then the yield y in Eq. (30) should be equal to the nucleosynthetic yield. However, galaxies do not evolve as closed systems, so we expect the overall metal yield to differ from the nucleosynthetic value. The effective yield y_{eff} can be defined by inverting Eq. (30) (e.g., Lequeux et al. 1979):
The effective yield is expected to be lower than the nucleosynthetic value y because metals in galaxies are lost from the system through outflows, and the metal content of the ISM gas will be diluted with infall of metalpoor gas (e.g., Edmunds 1990).
The effective yields y_{eff} for MAGMA, computed using Eq. (31), are shown in Fig. 7, plotted against baryonic gas mass fraction μ_{g}. The model prediction from Eq. (15) through our parameterizations of Δη and ζ_{w} is shown as two curves, corresponding to the two different Δη priors. For reference, our adopted value of y is shown as a horizontal shortdashed line. As also evident in previous figures, the results from the Δη ≥ 0 prior (longdashed orange curve) are in slightly better agreement with the data. The inflection of y_{eff} at high μ_{g} is due mainly to the strong downward curvature of Z_{g} as seen in Fig. 3. Moreover, y_{eff} is always (with a single possible exception at high μ_{g}, NGC 4731) significantly smaller than the true yield y = 0.037 adopted here. Edmunds (1990) has shown that this behavior for y_{eff} would be expected if metals were lost through stellar and SNedriven bulk outflows, or diluted by accretion of metalpoor gas. MAGMA shows a strong confirmation of this scenario. It is clear from Fig. 7 that our parameterized model is a good approximation to the data (see also Fig. 3).
Fig. 7. Effective yield as a function of μ_{g}. Galaxies are coded by gas fraction, F_{g}, as shown in the vertical color bar. MAGMA medians are shown as a heavy gray line, with ±1σ excursions as the surrounding gray regions. The two dashed curves correspond to the predictions from the Δη and ζ_{w} parameterizations (see also Fig. 3), with the symmetric Δη priors given as a shortdashed (purple) curve, and the Δη ≥ 0 priors as a longdashed (orange) one. The horizontal dashed line corresponds to our adopted true yield y = 0.037 (see text). 
Models that do not take into account the metalenriched bulk outflows at low M_{star}, namely assuming that Z_{w} = Z_{g}, predict that at high μ_{g}, y_{eff} would approach the nucleosynthetic value y (e.g., Erb 2008). This can be seen by expanding Eq. (31) in a Taylor’s series (see Dalcanton 2007; Filho et al. 2013) so that for high μ_{g}, we have:
Thus, in very gasrich galaxies, the effective yield would be approximately independent of gas fraction. As pointed out by Dalcanton (2007), metalpoor accreted gas will lower the metallicity of the system, but it will also augment μ_{g}, thus maintaining the system in a pseudoclosedbox equilibrium. The linear dependence of y_{eff} on metal mass M_{Z}, as shown by Eq. (32), also means that metalenriched outflows are very effective at lowering the effective yield in galaxies with high gasmass fractions (see Dalcanton 2007). Indeed, Fig. 7 shows that even at high μ_{g}, y_{eff} is significantly lower than the true yield y.
It is important to remember that Eq. (31) measures Z_{g} in the ionized gas, but the gas content μ_{g} is derived from neutral gas (H I + H_{2}). There is some evidence that the metal contents of these two gas phases are different; in UV absorption studies of metalpoor dwarf galaxies (e.g., Thuan et al. 2005; Lebouteiller et al. 2013), the metallicity of the neutral gas in these systems was always found to be lower than in the ionized gas. A similar conclusion was reached by Filho et al. (2013) and Thuan et al. (2016) by comparing y_{eff} with y, and attributing the difference to the lower metallicity in the neutral gas.
Figure 8 plots the effective yield y_{eff} against stellar mass M_{star} (left panel) and circular velocity, V_{circ} (right). Even though there is much scatter at low M_{star} and V_{circ} (M_{star} ≲ 3 × 10^{8} M_{⊙} and V_{circ} ≲ 50 km s^{−1}), y_{eff} is increasing with mass and velocity up to 3 × 10^{9} M_{⊙} and ∼100 km s^{−1}, in qualitative agreement with similar trends reported in the literature (e.g., Garnett 2002; Pilyugin et al. 2004; Tremonti et al. 2004; Dalcanton 2007; Ekta & Chengalur 2010). For more massive galaxies, there is an inflection in the relations with a mild decline at masses and velocities larger than 3 × 10^{9} M_{⊙} and ∼100 km s^{−1}.
Fig. 8. Effective yield as a function of M_{star} (left panel) and V_{circ} (right panel). Galaxies are coded by gas fraction, F_{g}, as shown in the vertical color bar. MAGMA medians are shown as a heavy gray line, with ±1σ excursions as the surrounding gray regions. Also shown in the left panel as two dashed curves are the predictions from the Δη and ζ_{w} parameterizations (see also Fig. 3), with the symmetric Δη priors given as a shortdashed (purple) curve, and the Δη ≥ 0 priors as a longdashed (orange) one. The right panel illustrates galaxies with high starformation efficiency (short gas depletion times) with a blue dashed curve, and those with low efficiency (long depletion times) as a red one. The horizontal shortdashed line corresponds to the yield of 0.037 adopted in this paper. See text for more details. 
The right panel of Fig. 8 also shows the MAGMA sample divided into two regimes in gas depletion times, τ_{g}, defined as M_{g}/SFR (see Paper II); galaxies with long τ_{g} (τ_{g} > 4 × 10^{9} yr) are shown as a dashed red curve, and short times τ_{g} < 4 × 10^{9} yr as a dashed blue one. The gas depletion time τ_{g} is the inverse of the socalled starformation efficiency (SFE, see Paper II), that also enters into Eq. (13), so short τ_{g} would correspond to higher SFE, and vice versa. There can be as much as a factor of ten in y_{eff} between these two categories of galaxies, with a mean difference of ∼2 − 3 in y_{eff}. At fixed mass, larger yields correspond to longer τ_{g} or lower SFE. In these galaxies with longer τ_{g}, the peak in y_{eff} at M_{star} ∼ 3 × 10^{9} M_{⊙} seems to disappear, while there is a steady decline of y_{eff} toward lower M_{star} for shorter depletion times (see also LaraLópez et al. 2019). Longer τ_{g} are also associated with higher gas fractions (see color coding in Fig. 8); thus Eq. (32) may be telling us that the higher effective yields are simply a function of the enhanced impact of metalenriched outflows in gasrich galaxies.
7. Characterizing stellar feedback through winds
By adopting the formalism for galactic chemical evolution formulated by Pagel (2009), we have quantified with MAGMA the metalweighted mass loading factor ζ_{w} (=η_{w}Z_{w}/Z_{g}), and with somewhat less success, the difference in mass loading of gas accretion and outflows, Δη. However, we are interested in the individual massloading factor relative to the SFR, η_{w}, and the metal content of the winds, Z_{w}. In this section, we explore the constraints for mass loading η_{w} and Z_{w} provided by our results.
7.1. The entrainment fraction
Assuming that Z_{w} = Z_{g}, as done in much previous work (e.g., Erb 2008; Finlator & Davé 2008; Dayal et al. 2013; Hunt et al. 2016; Belfiore et al. 2019), implies that virtually all of the material in the outflow is ISM gas. The gas metallicity could be only a lower limit to the metallicity of the wind (e.g., Peeples & Shankar 2011), because the metals produced by SNe would be expected to enrich the gas above that of the galaxy (e.g., Vader 1986). Thus, to estimate Z_{w}, we also need to consider the fraction of gas that has been entrained in the wind, f_{ent}, and the metallicity of the SNe ejecta, Z_{eject} (see Dalcanton 2007; Peeples & Shankar 2011):
In lowmass galaxies, f_{ent} is found to be low, both observationally (e.g., Martin et al. 2002) and from simulations (e.g., Mac Low & Ferrara 1999). The entrainment fraction is closely related to the massloading factor of the wind η_{w}, and probably varies with galaxy mass (gravitational potential, see e.g., Hopkins et al. 2012). The maximum oxygen abundance of SNe ejecta Z_{eject} is known to within a factor of 2−3, Z_{eject} ∼ 0.04 − 0.1 (e.g., Woosley & Weaver 1995), and can constrain f_{ent} and η_{w} (e.g., Martin et al. 2002). The increase of ζ_{w} with decreasing virial velocity and stellar mass (see Fig. 2) is possibly due to diminishing f_{ent} to raise Z_{w} closer to Z_{eject}. To verify this, we need to isolate η_{w}, the massloading factor in the product ζ_{w} = η_{w}Z_{w}/Z_{g}.
7.2. The driving mechanisms behind stellar outflows
One way to do this is to compare the V_{circ} scaling of ζ_{w} with expected V_{circ} scalings of η_{w} (see Eq. (22)). Galactic outflows powered by stellar feedback are thought to be driven either by momentum, imparted to the ISM by massive stellar winds and SNe through radiation pressure (e.g., Murray et al. 2005, 2010; Hopkins et al. 2011; Zhang & Thompson 2012), or by energy, injected into the ISM by massive stars and corecollapse SNe (e.g., Dekel & Silk 1986; Heckman et al. 1990; Dalla Vecchia & Schaye 2012; Christensen et al. 2016; Naab & Ostriker 2017).
For outflows that conserve momentum, we would expect a powerlaw scaling while for energyconserving winds, (e.g., Dekel & Silk 1986; Murray et al. 2005; Hopkins et al. 2012). The only way such scalings can be reconciled with the (modified powerlaw) form of Eq. (22) is by a complementary scaling of Z_{w}. Examples of such complementary scaling are shown in Fig. 9 where the ratio of wind and ISM metallicities Z_{w} and Z_{g} is plotted against virial velocity V_{circ}. In the left panel of Fig. 9, we have taken the normalization of Z_{w}/Z_{g} to be the same as that as the bestfit of Eq. (22) using our Bayesian estimates, simply V_{circ}/V_{0}, where V_{0} = 32 km s^{−1} for symmetric Δη priors, and V_{0} = 35 km s^{−1} for priors with Δη ≥ 0. Both fits have powerlaw index b ∼ 2. The assumptions of momentumconserving or energyconserving outflows lead to different estimates for the trend of Z_{w}/Z_{g}.
Fig. 9. Ratio of inferred wind and ISM metallicities Z_{w}/Z_{g} plotted against virial velocity V_{circ}. Also shown (with arbitrary scaling) is the shape of Bayesian ζ_{w} as also plotted in Fig. 2; it is clearly not a power law, thus the different powerlaw forms of η_{w} introduce inflections in Z_{w}/Z_{g}. Left panel: as indicated in the legend, different trends of Z_{w}/Z_{g} are inferred by dividing the fitted form of ζ_{w} (see Eq. (22)) by various scalings of η_{w} (see text for more details). The shaded regions correspond to the trends allowed by the two different Bayesian estimates of ζ_{w} according to the different Δη priors. Right panel: brokenpower law trend for η_{w} found by Muratov et al. (2015) is shown as a solid heavy (gray) curve, where the inflection at V_{circ} = 60 km s^{−1} corresponds to a sharp decrease of Z_{w}/Z_{g} (see text for more details). The shortdashed (darkblue) curve shows the relation found by Mitchell et al. (2020), in the lowmass regime where stellar feedback dominates, and the longdashed (blue) curve shows the observed relation derived by Chisholm et al. (2018). 
Because dwarf galaxies are found to have a low entrainment fraction (e.g., Mac Low & Ferrara 1999; Martin et al. 2002), we would expect Z_{w}/Z_{g} to increase toward low V_{circ} because low f_{ent} would enrich the outflows toward the higher superSolar abundance of the SNe ejecta, Z_{eject} (see Eq. (33)). Energyconserving winds result in a flattening of Z_{w}/Z_{g} toward low V_{circ} (blue curves), and a cubic trend of η_{w} with V_{circ} gives a strongly declining Z_{w}/Z_{g} (gray curves). The trend of is the only one that gives a rising trend of Z_{w}/Z_{g} at low V_{circ}. Thus, momentumconserving winds would be favored by our results.
Quantifying this conclusion is, however, dependent on the metallicity calibration. For the PP04N2 calibration used here, the roughly quadratic lowmass V_{circ} dependence is shallower than the cubic (or quartic) dependence found by Peeples & Shankar (2011) for different O/H calibrations. Other calibrations such as those by Tremonti et al. (2004) or Denicoló et al. (2002) give differing trends (with powerlaw slopes b > 3 at low V_{circ}), that are inconclusive for determining the nature of the massloading factor η_{w}. The PP04N2 calibration adopted here results in the shallowest increase of ζ_{w} toward low V_{circ} (powerlaw slope b ≈ 2), relative to other O/H calibrations, and this is reflected in a more incisive lowmass (low V_{circ}) increase in Z_{w}/Z_{g}.
The right panel of Fig. 9 shows the predictions for η_{w} from the subgrid FIRE (Feedback in Realistic Environments) simulations (Muratov et al. 2015), and the EAGLE simulations (Mitchell et al. 2020). At low V_{circ} (< 60 km s^{−1}), Muratov et al. (2015) find the powerlaw dependence on V_{circ} to be steep (∼ − 3.2), while at higher masses (V_{circ}) the winds are momentumdriven with η_{w} ∝ . A similar conclusion was reached by Davé et al. (2013), namely that at low masses, winds are energydriven (or even more steeply with V_{circ}), while momentumdriven winds dominate at higher masses. Judging from Fig. 9, these trends would be incompatible with our results. In contrast, EAGLE simulations (e.g., Mitchell et al. 2020) indicate a slope midway between, namely ∼ − 1.5, leading to a very gradual increase in metallicities toward lower V_{circ}. Comparing the left and righthand sides of Fig. 9 shows that such a dependence better reflects our results at low V_{circ}. At higher V_{circ} (masses), the situation is more complex because the flattening of ζ_{w} is difficult to interpret in terms of a single powerlaw dependence of η_{w} on V_{circ}. In any case, we conclude that at low masses, winds in stellar outflows cannot be energy driven, but rather are more consistent with momentumdriven outflows (), or a slightly steeper “combined” slope as found by the EAGLE models with .
There is also the question, not touched upon here, of whether or not outflows are able to expel material from the galaxy halo, thus enriching the IGM. The entrainment fraction f_{ent} would be key to relating the strength of the outflows with metallicity in the winds (e.g., Martin et al. 2002; Chisholm et al. 2018) and deserves further study. In addition, entrainment is also affected by the gas phases involved in the outflow. Our analysis traces the metallicity in the ionized gas (considered the “cool” phase by most simulations, T < 2 × 10^{4} K), thought to carry most of the mass, whereas physically realistic winds also contain a “hot phase” (T > 5 × 10^{5} K) that harbors most of the energy (Kim et al. 2020). Consequently, with MAGMA, we are missing the diagnostic of the hot Xray phase.
8. Summary and future perspectives
In this paper, we explore the impact of ejective feedback on galaxy evolution using data from the MAGMA sample of 392 galaxies with H I and CO detections, which collects a homogeneous sample of local field starforming galaxies, covering the widest range of stellar and gas masses, metalliticies and star formation compiled up to now. After presenting the sample, and the correlations among metallicity, SFR, M_{star} and (cool) gasphase components in Paper I and Paper II, in this third paper of the series, we investigate outflow metal loading factors and effective yields in terms of gas fraction, mass and circular velocity, using the formalism of Pagel (2009). The main results of the paper follow.

We model the chemical enrichment of MAGMA galaxies using a formalism based on Pagel (2009), without any constraint on the M_{g} time evolution and including the assumption that wind metallicity Z_{w} is different from the ISM metallicity Z_{g}. Wind metallicity Z_{w} can be greater than the ambient Z_{g} if it is primarily comprised of supernova ejecta or depressed if a sufficient amount of metalpoor gas is entrained as the wind propagates out of the galaxy (Chisholm et al. 2015; Creasey et al. 2015; Muratov et al. 2017; Pillepich et al. 2018). Using this formalism we have observationally constrained two parameters important for stellar feedback: (a) the metalloading factor of the winds ζ_{w}; and (b) the difference of the mass loading of the accretion and winds, Δη = η_{a} − η_{w}. This is done through a Bayesian approach applied to subsamples of galaxies in different M_{star} and μ_{g} bins.

In terms of trends with M_{star}, we find that Δη is degenerate, and exceedingly difficult to constrain. However, our Bayesian estimates based on μ_{g} bins show that Δη increases significantly with increasing gas fraction, implying that Δη is a strong function of gas content and possibly evolutionary state (see Appendix B).

The wind metalloading factors ζ_{w} are found to increase with decreasing stellar mass and circular velocity, with a flattening toward larger M_{star}. This indicates that the efficiency of metalloading in the winds depends on the potential depth and is strongly enhanced in lowmass systems, in agreement with previous work based on different methodologies for estimating ζ_{w} (e.g., Peeples & Shankar 2011; Chisholm et al. 2018).

We have also constrained effective yields y_{eff} with MAGMA data, and find that y_{eff} increases with stellar mass and V_{circ} up to a certain threshold, followed by a flattening toward larger M_{star}, and possibly a mild inversion at the highest masses.

Finally, we have inferred possible trends with V_{circ} for the wind massloading factors η_{w} and find that our results for low V_{circ} (M_{star}) favor momentumdriven winds, rather than energydriven ones.
Various quantities in our analysis, ζ_{w}, y_{eff}, and Z_{g}, present inflections at M_{star} ∼ 3 × 10^{9} M_{⊙}, corresponding to virial velocities ∼100 km s^{−1}. This corresponds to a “gasrichness” mass scale, associated with transitions of physical processes, separating two different mass regimes (Paper II; see also Kannappan 2004). It delineates an accretiondominated regime at low masses, where stellar outflows are relevant, and an intermediate mass, gasequilibrium regime, where outflows are less effective, and galaxies are able to consume gas through SF as fast as it is accreted. There would also be an additional, larger “bimodality” mass threshold at M_{star} ∼ 3 × 10^{10} M_{⊙} (Paper II, Kannappan et al. 2013) but this is only weakly evident in MAGMA, mainly because of the exclusion of passive or “quenched” gaspoor galaxies and AGN.
These results can also be interpreted within a broader context, since the mass thresholds defining metallicity and gas properties are associated with changes of correlations involving stellar populations and stellar and dark matter distributions (Tortora et al. 2010, 2019). In particular, inflections with mass of metalloading and metal yields at the “gasrichness” threshold seem to define changes of the slopes of total mass density profiles including dark matter (e.g., Tortora et al. 2019); this threshold approximatively separates galaxies with mass density slopes shallower and steeper than ∼ − 1 (see their Fig. 1). More investigations on this mass scale will follow.
As outlined in Paper II, further improvement of the MAGMA coverage, with more detections of total gas content in dwarf galaxies, will strengthen the conclusions of our series of papers. In addition, it is important to reassess the impact on our conclusions of assumptions made here, namely that the accretion rate Ṁ_{a} is proportional to the SFR and the chemically unenriched nature of the accreted gas. In this paper, we have also pursued a first step to combine information from completely independent datasets and observables. Only by including gas, stellar, kinematical and dynamical observables, will we be able to determine the efficiency of the physical phenomena which shape and regulate galaxy evolution.
For all statistical analysis, we rely on the R statistical package: R Core Team (2018), R: A language and environment for statistical computing, R Foundation for Statistical Computing, Vienna, Austria (https://www.Rproject.org/).
Acknowledgments
We thank the anonymous referee and the editor for profound insight that helped to greatly improve the paper. We are grateful for funding from the INAF PRINSKA 2017 program 1.05.01.88.04.
References
 Andrews, B. H., & Martini, P. 2013, ApJ, 765, 140 [NASA ADS] [CrossRef] [Google Scholar]
 Baldwin, J. A., Phillips, M. M., & Terlevich, R. 1981, PASP, 93, 5 [Google Scholar]
 Behroozi, P. S., Wechsler, R. H., & Conroy, C. 2013, ApJ, 770, 57 [NASA ADS] [CrossRef] [Google Scholar]
 Belfiore, F., Vincenzo, F., Maiolino, R., & Matteucci, F. 2019, MNRAS, 487, 456 [NASA ADS] [CrossRef] [Google Scholar]
 Boselli, A., Boissier, S., Cortese, L., et al. 2009, ApJ, 706, 1527 [Google Scholar]
 Boselli, A., Cortese, L., & Boquien, M. 2014a, A&A, 564, A65 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Boselli, A., Cortese, L., Boquien, M., et al. 2014b, A&A, 564, A66 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Brinchmann, J., Charlot, S., White, S. D. M., et al. 2004, MNRAS, 351, 1151 [Google Scholar]
 Catinella, B., Saintonge, A., Janowiecki, S., et al. 2018, MNRAS, 476, 875 [NASA ADS] [CrossRef] [Google Scholar]
 Chabrier, G. 2003, PASP, 115, 763 [Google Scholar]
 Chisholm, J., Tremonti, C. A., Leitherer, C., et al. 2015, ApJ, 811, 149 [NASA ADS] [CrossRef] [Google Scholar]
 Chisholm, J., Tremonti, C., & Leitherer, C. 2018, MNRAS, 481, 1690 [NASA ADS] [CrossRef] [Google Scholar]
 Christensen, C. R., Davé, R., Governato, F., et al. 2016, ApJ, 824, 57 [CrossRef] [Google Scholar]
 Creasey, P., Theuns, T., & Bower, R. G. 2015, MNRAS, 446, 2125 [CrossRef] [Google Scholar]
 Dalcanton, J. J. 2007, ApJ, 658, 941 [CrossRef] [Google Scholar]
 Dalla Vecchia, C., & Schaye, J. 2012, MNRAS, 426, 140 [NASA ADS] [CrossRef] [Google Scholar]
 Davé, R., Finlator, K., & Oppenheimer, B. D. 2012, MNRAS, 421, 98 [Google Scholar]
 Davé, R., Katz, N., Oppenheimer, B. D., Kollmeier, J. A., & Weinberg, D. H. 2013, MNRAS, 434, 2645 [CrossRef] [Google Scholar]
 Dayal, P., Ferrara, A., & Dunlop, J. S. 2013, MNRAS, 430, 2891 [Google Scholar]
 Dekel, A., & Birnboim, Y. 2006, MNRAS, 368, 2 [NASA ADS] [CrossRef] [Google Scholar]
 Dekel, A., & Silk, J. 1986, ApJ, 303, 39 [Google Scholar]
 Denicoló, G., Terlevich, R., & Terlevich, E. 2002, MNRAS, 330, 69 [CrossRef] [Google Scholar]
 De Rossi, M. E., Bower, R. G., Font, A. S., Schaye, J., & Theuns, T. 2017, MNRAS, 472, 3354 [Google Scholar]
 Edmunds, M. G. 1990, MNRAS, 246, 678 [NASA ADS] [Google Scholar]
 Edmunds, M. G., & Pagel, B. E. J. 1984, in Three Problems in the Chemical Evolution of Galaxies, eds. C. Chiosi, & A. Renzini, Astrophys. Space Sci. Lib., 109, 341 [NASA ADS] [Google Scholar]
 Ekta, B., & Chengalur, J. N. 2010, MNRAS, 406, 1238 [NASA ADS] [Google Scholar]
 Elmegreen, B. G., Rubio, M., Hunter, D. A., et al. 2013, Nature, 495, 487 [NASA ADS] [CrossRef] [Google Scholar]
 Emami, N., Siana, B., Weisz, D. R., et al. 2019, ApJ, 881, 71 [NASA ADS] [CrossRef] [Google Scholar]
 Erb, D. K. 2008, ApJ, 674, 151 [NASA ADS] [CrossRef] [Google Scholar]
 Erb, D. K., Shapley, A. E., Pettini, M., et al. 2006, ApJ, 644, 813 [NASA ADS] [CrossRef] [Google Scholar]
 Filho, M. E., Winkel, B., Sánchez Almeida, J., et al. 2013, A&A, 558, A18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Finlator, K., & Davé, R. 2008, MNRAS, 385, 2181 [NASA ADS] [CrossRef] [Google Scholar]
 Fraternali, F. 2017, in Gas Accretion via Condensation and Fountains, eds. A. Fox, & R. Davé, Astrophys. Space Sci. Lib., 430, 323 [Google Scholar]
 Fraternali, F., & Tomassetti, M. 2012, MNRAS, 426, 2166 [NASA ADS] [CrossRef] [Google Scholar]
 Garnett, D. R. 2002, ApJ, 581, 1019 [NASA ADS] [CrossRef] [Google Scholar]
 Ginolfi, M., Hunt, L. K., Tortora, C., Schneider, R., & Cresci, G. 2020, A&A, 638, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Graziani, L., Salvadori, S., Schneider, R., et al. 2015, MNRAS, 449, 3137 [NASA ADS] [CrossRef] [Google Scholar]
 Graziani, L., de Bennassuti, M., Schneider, R., Kawata, D., & Salvadori, S. 2017, MNRAS, 469, 1101 [NASA ADS] [CrossRef] [Google Scholar]
 Haynes, M. P., Giovanelli, R., Martin, A. M., et al. 2011, AJ, 142, 170 [Google Scholar]
 Heckman, T. M., Armus, L., & Miley, G. K. 1990, ApJS, 74, 833 [Google Scholar]
 Hopkins, P. F., Quataert, E., & Murray, N. 2011, MNRAS, 417, 950 [NASA ADS] [CrossRef] [Google Scholar]
 Hopkins, P. F., Quataert, E., & Murray, N. 2012, MNRAS, 421, 3522 [Google Scholar]
 Huang, S., Haynes, M. P., Giovanelli, R., & Brinchmann, J. 2012, ApJ, 756, 113 [Google Scholar]
 Hunt, L., Magrini, L., Galli, D., et al. 2012, MNRAS, 427, 906 [NASA ADS] [CrossRef] [Google Scholar]
 Hunt, L., Dayal, P., Magrini, L., & Ferrara, A. 2016, MNRAS, 463, 2002 [NASA ADS] [CrossRef] [Google Scholar]
 Hunt, L. K., De Looze, I., Boquien, M., et al. 2019, A&A, 621, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hunt, L. K., Tortora, C., Ginolfi, M., & Schneider, R. 2020, A&A, 643, A180 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Izotov, Y. I., Thuan, T. X., & Stasińska, G. 2007, ApJ, 662, 15 [NASA ADS] [CrossRef] [Google Scholar]
 Kannappan, S. J. 2004, ApJ, 611, L89 [Google Scholar]
 Kannappan, S. J., Stark, D. V., Eckert, K. D., et al. 2013, ApJ, 777, 42 [Google Scholar]
 Kewley, L. J., & Ellison, S. L. 2008, ApJ, 681, 1183 [Google Scholar]
 Kim, C.G., Ostriker, E. C., Somerville, R. S., et al. 2020, ApJ, 900, 61 [Google Scholar]
 Krumholz, M. R., Fumagalli, M., da Silva, R. L., Rendahl, T., & Parra, J. 2015, MNRAS, 452, 1447 [NASA ADS] [CrossRef] [Google Scholar]
 LaraLópez, M. A., Cepa, J., Bongiovanni, A., et al. 2010, A&A, 521, L53 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 LaraLópez, M. A., De Rossi, M. E., Pilyugin, L. S., et al. 2019, MNRAS, 490, 868 [CrossRef] [Google Scholar]
 Larson, R. B. 1972, Nat. Phys. Sci., 236, 7 [NASA ADS] [CrossRef] [Google Scholar]
 Lebouteiller, V., Heap, S., Hubeny, I., & Kunth, D. 2013, A&A, 553, A16 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lelli, F., McGaugh, S. S., & Schombert, J. M. 2016a, ApJ, 816, L14 [Google Scholar]
 Lelli, F., McGaugh, S. S., & Schombert, J. M. 2016b, AJ, 152, 157 [Google Scholar]
 Lequeux, J., Peimbert, M., Rayo, J. F., Serrano, A., & TorresPeimbert, S. 1979, A&A, 500, 145 [NASA ADS] [Google Scholar]
 Leroy, A. K., Walter, F., Brinks, E., et al. 2008, AJ, 136, 2782 [Google Scholar]
 Leroy, A. K., Walter, F., Bigiel, F., et al. 2009, AJ, 137, 4670 [Google Scholar]
 Leroy, A. K., Sandstrom, K. M., Lang, D., et al. 2019, ApJS, 244, 24 [Google Scholar]
 Lilly, S. J., Carollo, C. M., Pipino, A., Renzini, A., & Peng, Y. 2013, ApJ, 772, 119 [NASA ADS] [CrossRef] [Google Scholar]
 Ma, X., Hopkins, P. F., FaucherGiguère, C.A., et al. 2016, MNRAS, 456, 2140 [NASA ADS] [CrossRef] [Google Scholar]
 Mac Low, M.M., & Ferrara, A. 1999, ApJ, 513, 142 [NASA ADS] [CrossRef] [Google Scholar]
 Madau, P., Weisz, D. R., & Conroy, C. 2014, ApJ, 790, L17 [NASA ADS] [CrossRef] [Google Scholar]
 Maiolino, R., Nagao, T., Grazian, A., et al. 2008, A&A, 488, 463 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mannucci, F., Cresci, G., Maiolino, R., Marconi, A., & Gnerucci, A. 2010, MNRAS, 408, 2115 [NASA ADS] [CrossRef] [Google Scholar]
 Martin, C. L., Kobulnicky, H. A., & Heckman, T. M. 2002, ApJ, 574, 663 [NASA ADS] [CrossRef] [Google Scholar]
 Mitchell, P. D., Schaye, J., Bower, R. G., & Crain, R. A. 2020, MNRAS, 494, 3971 [Google Scholar]
 Mitra, S., Davé, R., & Finlator, K. 2015, MNRAS, 452, 1184 [NASA ADS] [CrossRef] [Google Scholar]
 Moster, B. P., Somerville, R. S., Maulbetsch, C., et al. 2010, ApJ, 710, 903 [Google Scholar]
 Moster, B. P., Naab, T., & White, S. D. M. 2013, MNRAS, 428, 3121 [Google Scholar]
 Muratov, A. L., Kereš, D., FaucherGiguère, C.A., et al. 2015, MNRAS, 454, 2691 [NASA ADS] [CrossRef] [Google Scholar]
 Muratov, A. L., Kereš, D., FaucherGiguère, C.A., et al. 2017, MNRAS, 468, 4170 [NASA ADS] [CrossRef] [Google Scholar]
 Murray, N., Quataert, E., & Thompson, T. A. 2005, ApJ, 618, 569 [NASA ADS] [CrossRef] [Google Scholar]
 Murray, N., Quataert, E., & Thompson, T. A. 2010, ApJ, 709, 191 [NASA ADS] [CrossRef] [Google Scholar]
 Naab, T., & Ostriker, J. P. 2017, ARA&A, 55, 59 [Google Scholar]
 Pagel, B. E. J. 2009, Nucleosynthesis and Chemical Evolution of Galaxies (Cambridge: Cambridge University Press) [CrossRef] [Google Scholar]
 Pagel, B. E. J., & Edmunds, M. G. 1981, ARA&A, 19, 77 [NASA ADS] [CrossRef] [Google Scholar]
 Pagel, B. E. J., & Patchett, B. E. 1975, MNRAS, 172, 13 [NASA ADS] [CrossRef] [Google Scholar]
 Peeples, M. S., & Shankar, F. 2011, MNRAS, 417, 2962 [NASA ADS] [CrossRef] [Google Scholar]
 Peng, Y.J., & Maiolino, R. 2014, MNRAS, 438, 262 [NASA ADS] [CrossRef] [Google Scholar]
 Pettini, M., & Pagel, B. E. J. 2004, MNRAS, 348, L59 [NASA ADS] [CrossRef] [Google Scholar]
 Pillepich, A., Springel, V., Nelson, D., et al. 2018, MNRAS, 473, 4077 [Google Scholar]
 Pilyugin, L. S., Contini, T., & Vílchez, J. M. 2004, A&A, 423, 427 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Recchi, S., Spitoni, E., Matteucci, F., & Lanfranchi, G. A. 2008, A&A, 489, 555 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Romano, D., Karakas, A. I., Tosi, M., & Matteucci, F. 2010, A&A, 522, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rubio, M., Elmegreen, B. G., Hunter, D. A., et al. 2015, Nature, 525, 218 [NASA ADS] [CrossRef] [Google Scholar]
 Saintonge, A., Kauffmann, G., Kramer, C., et al. 2011, MNRAS, 415, 32 [NASA ADS] [CrossRef] [Google Scholar]
 Saintonge, A., Catinella, B., Tacconi, L. J., et al. 2017, ApJS, 233, 22 [Google Scholar]
 Sales, L. V., Navarro, J. F., Oman, K., et al. 2017, MNRAS, 464, 2419 [NASA ADS] [CrossRef] [Google Scholar]
 Sánchez Almeida, J., Elmegreen, B. G., MuñozTuñón, C., & Elmegreen, D. M. 2014, A&ARv, 22, 71 [Google Scholar]
 Schruba, A., Leroy, A. K., Walter, F., et al. 2012, AJ, 143, 138 [NASA ADS] [CrossRef] [Google Scholar]
 Searle, L., & Sargent, W. L. W. 1972, ApJ, 173, 25 [NASA ADS] [CrossRef] [Google Scholar]
 Shi, Y., Wang, J., Zhang, Z.Y., et al. 2015, ApJ, 804, L11 [CrossRef] [Google Scholar]
 Shi, Y., Wang, J., Zhang, Z.Y., et al. 2016, Nat. Commun., 7, 13789 [NASA ADS] [CrossRef] [Google Scholar]
 Spitoni, E., Calura, F., Matteucci, F., & Recchi, S. 2010, A&A, 514, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Tacchella, S., Dekel, A., Carollo, C. M., et al. 2016, MNRAS, 457, 2790 [Google Scholar]
 Thuan, T. X., Lecavelier des Etangs, A., & Izotov, Y. I. 2005, ApJ, 621, 269 [NASA ADS] [CrossRef] [Google Scholar]
 Thuan, T. X., Goehring, K. M., Hibbard, J. E., Izotov, Y. I., & Hunt, L. K. 2016, MNRAS, 463, 4268 [NASA ADS] [CrossRef] [Google Scholar]
 Tinsley, B. M. 1980, Fund. Cosm. Phys., 5, 287 [Google Scholar]
 Torrey, P., Vogelsberger, M., Marinacci, F., et al. 2019, MNRAS, 484, 5587 [NASA ADS] [Google Scholar]
 Tortora, C., Napolitano, N. R., Cardone, V. F., et al. 2010, MNRAS, 407, 144 [Google Scholar]
 Tortora, C., Posti, L., Koopmans, L. V. E., & Napolitano, N. R. 2019, MNRAS, 489, 5483 [Google Scholar]
 Tosi, M. 1982, ApJ, 254, 699 [NASA ADS] [CrossRef] [Google Scholar]
 Tremonti, C. A., Heckman, T. M., Kauffmann, G., et al. 2004, ApJ, 613, 898 [Google Scholar]
 Vader, J. P. 1986, ApJ, 305, 669 [NASA ADS] [CrossRef] [Google Scholar]
 Vincenzo, F., Matteucci, F., Belfiore, F., & Maiolino, R. 2016, MNRAS, 455, 4183 [Google Scholar]
 Weisz, D. R., Johnson, B. D., Johnson, L. C., et al. 2012, ApJ, 744, 44 [NASA ADS] [CrossRef] [Google Scholar]
 Woosley, S. E., & Weaver, T. A. 1995, ApJS, 101, 181 [Google Scholar]
 Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868 [Google Scholar]
 Wu, P.F., Zahid, H. J., Hwang, H. S., & Geller, M. J. 2017, MNRAS, 468, 1881 [NASA ADS] [CrossRef] [Google Scholar]
 Zahid, H. J., Dima, G. I., Kudritzki, R.P., et al. 2014, ApJ, 791, 130 [NASA ADS] [CrossRef] [Google Scholar]
 Zhang, D., & Thompson, T. A. 2012, MNRAS, 424, 1170 [NASA ADS] [CrossRef] [Google Scholar]
Appendix A: Additional cornerplots for Bayesian best fits
Here we show the additional corner plots of the χ^{2} surface as a function of the model priors (Δη, ζ_{w}) for MAGMA, for Δη ∈ [ − 10, 3.4] (Fig. A.1) and Δη ∈ [0, 3.4] (Fig. A.2).
Fig. A.1. Corner plots of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for MAGMA. The violet contours correspond to the minimum χ^{2} value. The top and right panels of each corner plot report the probability density distributions for the marginalized parameters; confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF mode for Δη, median for ζ_{w}) is shown by a vertical dashed line. Here, Δη ∈ [ − 10, 3.4] (see Fig. 1 for the “symmetric” Δη intervals). 
Fig. A.2. Corner plots of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for MAGMA. The violet contours correspond to the minimum χ^{2} value. The top and right panels of each corner plot report the probability density distributions for the marginalized parameters; confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF median) is shown by a vertical dashed line. Here, Δη ∈ [0,3.4] (see Fig. 1 for the “symmetric” Δη intervals). 
Appendix B: SFHs and toymodel predictions for gas and stellar growth
In order to investigate the time derivatives of M_{g} and M_{star} and their ratio, we introduce an arbitrary SFH. As a simple example, we have taken an exponential with arbitrary scaling defined by β and the starformation time scale τ_{SF}:
Equation (B.1) can be integrated to obtain M_{star} [see Eq. (2)], at the current age of the galaxy t_{gal}:
If we take the timescale of star formation τ_{SF} to be equal to the gas depletion time τ_{g} [the inverse of ϵ_{⋆}, see Eq. (13)], there are two unknowns in Eq. (B.2): t_{gal} and ψ_{0} (β is fixed arbitrarily for rising and declining SFHs). However, ψ(t_{gal}) = ψ_{0} exp(βt_{gal}/τ_{SF}) = SFR and M_{⋆}(t_{gal}) = M_{star}, are the observables that we measure, so that Eq. (B.2) can be solved for ψ_{0}:
This toy model also gives the initial gas mass M_{i} = ψ_{0}τ_{SF}.
The remaining variable needed to define the galaxy’s evolutionary state is t_{gal}, which can be established with ψ_{0} through Eq. (B.1):
In case of a constant SFH (β = 0), with ψ(t) = ψ_{0} = SFR, the solutions are quite simple: M_{⋆}(t) = M_{star}, M_{gas} = M_{i} = SFR τ_{SF}, and Δη = α.
Figure B.1 shows the values inferred for t_{gal} with the simple SFH of Eq. (B.1) with β = −1 (declining exponential) plotted against M_{star} (top panel) and ψ_{0} plotted against M_{star} (middle). Although there is significant scatter, galaxies with M_{star} ≲ 10^{9} M_{⊙} are apparently younger, with ages as low as 10^{9} yr, rather than the overall inferred age of more massive galaxies, 10^{10} yr. The initial SFR ψ_{0} is very tightly linked with M_{star}, over a dynamic range of ∼10^{3}. The bottom panel of Fig. B.1 illustrates that the total gas depletion times (i.e., including both molecular and atomic components) are relatively constant with M_{star}.
Fig. B.1. Inferred galaxy age t_{gal} plotted against M_{star} for MAGMA galaxies (top panel); initial SFR at t = 0, ψ_{0} plotted against M_{star} (middle); τ_{g} (assumed equal to τ_{SF}) plotted against M_{star} (bottom). In the two top panels, we assumed a simple declining exponential SFH with β = −1 [Eq. (B.1)], and that . 
Now the only missing ingredient to describe M_{g} is Δη [see Eq. (14)]. Ultimately, the toy model presented here requires through Eq. (13) that the SFH given by Eq. (B.1) is related to the time evolution of M_{g} [Eq. (14)]. If is the same as ϵ_{⋆}, as we have assumed, the formalism requires:
Thus, fixing β allows only certain values of Δη. Figure B.2 shows the time evolution of M_{g} and M_{star} over the last 50% of the galaxies’ lifetime, assuming that t_{gal} is given by Eq. (B.4), and ψ_{0} by Eq. (B.3). The different sets of curves correspond to different M_{star} medians and the corresponding medians of τ_{SF}, t_{gal}, and ψ_{0}, while the vertical separation is achieved through the assignment of different Δη as required by β = ( − 1, 0, 0.1). The rising exponential SFH cannot have β arbitrarily large because otherwise Eq. (B.3) has negative values; β = 0.1 is close to the largest value allowed by the data. Figure B.2 clearly shows that according to the SFH, and the value of β and consequently Δη, dlog M_{g}/dlog M_{star} can be positive or negative. Thus, it is difficult to assign Δη through the application of Eq. (27).
Fig. B.2. Time evolution of M_{g} vs M_{star} over the last 50% of the galaxies’ lifetime, according to the toymodel based on three exponential SFHs: β = ( − 1, 0, 1), and also assuming that the SFR timescale τ_{SF} is the same as the gas depletion time, τ_{g} = . Galaxy age t_{gal} is given by Eq. (B.4), ψ_{0} by Eq. (B.3), and the M_{g} evolution by Eq. (14). The median curves in different M_{star} bins are plotted, where the (median observed) gas mass M_{g} is taken to be M_{g}(t_{gal}). For Δη = α, the equilibrium solution, M_{g} = M_{i} throughout the lifetime of the galaxy. 
Appendix C: The formalism for the “homogeneous wind”
In this paper we have determined solutions for Z_{g} in the most general case of Z_{w} ≠ Z_{g}. We demonstrate (as already shown by Peeples & Shankar 2011) that in general this is a valid assumption since the outflowing gas Z_{w} is not the same as the ISM gas, Z_{g}, especially at low masses (e.g., Chisholm et al. 2015; Creasey et al. 2015; Muratov et al. 2017).
However, much earlier work assumed, instead, that Z_{w} = Z_{g}, adopting the socalled homogeneous wind model (e.g., Pagel 2009; Erb 2008). This assumption leads to ζ_{w} = η_{w}, and modifies Eq. (15) as follows:
This solution, already reported by Erb (2008) in her Eq. 11, places the dependence of η_{w} only in the exponent of M_{g}/M_{i}, while η_{a} also appears in the denominator in the Z_{g} expression.
The very interesting outcome of this exercise is that the η_{w} derived in this way reproduce the values of ζ_{w} shown in Fig. 3 and obtained using the more general formalism (using Z_{w} ≠ Z_{g}) of Eq. (15). Thus, metal loading ζ_{w} and mass loading η_{w} in the winds can be confused, depending on the assumptions made in the application of the formalism.
All Tables
All Figures
Fig. 1. Corner plots of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for MAGMA binned by M_{star}. The violet contours correspond to the minimum χ^{2} value. The top and right panels of each corner plot report the probability density distributions for the marginalized parameters; confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF median) is shown by a vertical dashed line. 

In the text 
Fig. 2. Difference of massloading factors Δη (top panel) and metalloading factors ζ_{w} (bottom panel) as a function of virial velocities obtained from the BTFR from Sales et al. (2017). In both panels, filled symbols show the Bayesian results, and in the top panel, open ones illustrate the upper bound of the PDF 1σ confidence intervals (CI). Different symbol types correspond to the interval of Δη in the Bayesian priors, as described in the text; the symbols for Δη ∈ [0, 3.4] and Δη ∈ [ − 10, 3.4] are (arbitrarily) offset along the abscissa for better visibility. The gray regions give the ±1σ uncertainties in the bestfit values. Also shown is the “equilibrium” asymptote, namely Δη = α (see Sect. 3.1). In the top panel, the two estimates of Δη from different priors are exactly coincident. In the bottom panel, the formulations by Peeples & Shankar (2011) for two O/H calibrations, redone here for our value of y, are given by (purple) shortdashed and (blue) longdashed curves for the PP04N2 and D02 calibrations, respectively. Also shown are the metalloading factors computed by Chisholm et al. (2018) from COS data. The text gives more details. 

In the text 
Fig. 3. Difference of massloading factors Δη and metalloading factors ζ_{w} inferred from the Bayesian analysis, and observed Z_{g} in the bottom panel, as a function of M_{star} on the left, and as a function of μ_{g} on the right. In the two upper panels, the bestfit Bayesian values are shown as filled symbols (circles or squares according to the Δη prior), with error bars corresponding to the ±1σ excursion, which is also illustrated by lightgray shaded regions. The rightmost point (magenta circle) in the top and middle right panels corresponds to the “extra” highμ_{g} bin shown in Fig. 6. In the top panels, a horizontal dashed line illustrates the “equilibrium solution” with Δη = α. In the top left panel, the (purple) dashed curve is a linear powerlaw fit of Δη to log(M_{star}), while the dashed curve in the top right panel shows Eq. (27), with d log M_{g}/d log M_{star} = 0.72, as found in Paper II. The fits to the curves for ζ_{w} in the middle panels are cubic: log(ζ_{w}) as a function of log(M_{star}) in the left panel, and ζ_{w} as a function of F_{g} in the right. In the lower panels, observed MAGMA Z_{g} are shown as small filled circles, color coded by gas fraction, while the Z_{g} medians and 1σ deviations are shown by large (black) filled circles and error bars. The dashed curves in the lower panels correspond to the calculation of Z_{g} based on the fits in the upper panels (as a function of μ_{g} on the left, and log(M_{star}) on the right): the longdashed curves (dark orange) correspond to the Δη ≥ 0 priors, while the shortdashed ones (purple) to the “symmetric” priors. The green curves in the left panels are described in the text, obtained by fitting individual galaxy Δη and ζ_{w} calculated using Eqs. (27) and (28). 

In the text 
Fig. 4. Baryonic gas fraction μ_{g} plotted against (log) M_{star} for MAGMA galaxies. Galaxies are coded by gas fraction, F_{g}, as shown in the vertical color bar. MAGMA medians are shown as a heavy gray line, with ±1σ excursions as the gray region. As noted in the text, the two lowmass galaxies with low μ_{g} are Sextans A and WLM, with probable substantial underestimates of the total gas content. 

In the text 
Fig. 5. Corner plots of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for MAGMA binned by μ_{g}. The violet contours correspond to the minimum χ^{2} value. The top and right panels of each corner plot report the probability density distributions for the marginalized parameters; confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF mode for Δη, median for ζ_{w}) is shown by a vertical dashed line. 

In the text 
Fig. 6. Corner plot of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for the highest μ_{g} bin MAGMA, here with a higher boundary than in Fig. 5: μ_{g} > 0.75. As in previous figures, the violet contours correspond to the minimum χ^{2} value, and the top and right panels report the probability density distributions for the marginalized parameters. Confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF mode for Δη, median for ζ_{w}) is shown by a vertical dashed line. 

In the text 
Fig. 7. Effective yield as a function of μ_{g}. Galaxies are coded by gas fraction, F_{g}, as shown in the vertical color bar. MAGMA medians are shown as a heavy gray line, with ±1σ excursions as the surrounding gray regions. The two dashed curves correspond to the predictions from the Δη and ζ_{w} parameterizations (see also Fig. 3), with the symmetric Δη priors given as a shortdashed (purple) curve, and the Δη ≥ 0 priors as a longdashed (orange) one. The horizontal dashed line corresponds to our adopted true yield y = 0.037 (see text). 

In the text 
Fig. 8. Effective yield as a function of M_{star} (left panel) and V_{circ} (right panel). Galaxies are coded by gas fraction, F_{g}, as shown in the vertical color bar. MAGMA medians are shown as a heavy gray line, with ±1σ excursions as the surrounding gray regions. Also shown in the left panel as two dashed curves are the predictions from the Δη and ζ_{w} parameterizations (see also Fig. 3), with the symmetric Δη priors given as a shortdashed (purple) curve, and the Δη ≥ 0 priors as a longdashed (orange) one. The right panel illustrates galaxies with high starformation efficiency (short gas depletion times) with a blue dashed curve, and those with low efficiency (long depletion times) as a red one. The horizontal shortdashed line corresponds to the yield of 0.037 adopted in this paper. See text for more details. 

In the text 
Fig. 9. Ratio of inferred wind and ISM metallicities Z_{w}/Z_{g} plotted against virial velocity V_{circ}. Also shown (with arbitrary scaling) is the shape of Bayesian ζ_{w} as also plotted in Fig. 2; it is clearly not a power law, thus the different powerlaw forms of η_{w} introduce inflections in Z_{w}/Z_{g}. Left panel: as indicated in the legend, different trends of Z_{w}/Z_{g} are inferred by dividing the fitted form of ζ_{w} (see Eq. (22)) by various scalings of η_{w} (see text for more details). The shaded regions correspond to the trends allowed by the two different Bayesian estimates of ζ_{w} according to the different Δη priors. Right panel: brokenpower law trend for η_{w} found by Muratov et al. (2015) is shown as a solid heavy (gray) curve, where the inflection at V_{circ} = 60 km s^{−1} corresponds to a sharp decrease of Z_{w}/Z_{g} (see text for more details). The shortdashed (darkblue) curve shows the relation found by Mitchell et al. (2020), in the lowmass regime where stellar feedback dominates, and the longdashed (blue) curve shows the observed relation derived by Chisholm et al. (2018). 

In the text 
Fig. A.1. Corner plots of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for MAGMA. The violet contours correspond to the minimum χ^{2} value. The top and right panels of each corner plot report the probability density distributions for the marginalized parameters; confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF mode for Δη, median for ζ_{w}) is shown by a vertical dashed line. Here, Δη ∈ [ − 10, 3.4] (see Fig. 1 for the “symmetric” Δη intervals). 

In the text 
Fig. A.2. Corner plots of χ^{2} surface as a function of the model parameters (Δη, ζ_{w}) for MAGMA. The violet contours correspond to the minimum χ^{2} value. The top and right panels of each corner plot report the probability density distributions for the marginalized parameters; confidence intervals (±1σ) are shown as violettinted shaded rectangular regions, and the MLE (PDF median) is shown by a vertical dashed line. Here, Δη ∈ [0,3.4] (see Fig. 1 for the “symmetric” Δη intervals). 

In the text 
Fig. B.1. Inferred galaxy age t_{gal} plotted against M_{star} for MAGMA galaxies (top panel); initial SFR at t = 0, ψ_{0} plotted against M_{star} (middle); τ_{g} (assumed equal to τ_{SF}) plotted against M_{star} (bottom). In the two top panels, we assumed a simple declining exponential SFH with β = −1 [Eq. (B.1)], and that . 

In the text 
Fig. B.2. Time evolution of M_{g} vs M_{star} over the last 50% of the galaxies’ lifetime, according to the toymodel based on three exponential SFHs: β = ( − 1, 0, 1), and also assuming that the SFR timescale τ_{SF} is the same as the gas depletion time, τ_{g} = . Galaxy age t_{gal} is given by Eq. (B.4), ψ_{0} by Eq. (B.3), and the M_{g} evolution by Eq. (14). The median curves in different M_{star} bins are plotted, where the (median observed) gas mass M_{g} is taken to be M_{g}(t_{gal}). For Δη = α, the equilibrium solution, M_{g} = M_{i} throughout the lifetime of the galaxy. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.