SDSS-IV MaNGA: drivers of stellar metallicity in nearby galaxies

ABSTRACT

The distribution of stellar metallicities within and across galaxies is an excellent relic of the chemical evolution across cosmic time. We present a detailed analysis of spatially resolved stellar populations based on >2.6 million spatial bins from 7439 nearby galaxies in the Sloan Digital Sky Survey-IV (SDSS-IV) Mapping Nearby Galaxies at APO (MaNGA) survey. To account for accurate inclination corrections, we derive an equation for morphology-dependent determination of galaxy inclinations. Our study goes beyond the well-known global mass–metallicity relation and radial metallicity gradients by providing a statistically sound exploration of local relations between stellar metallicity [Z/H], stellar surface mass density Σ_⋆, and galactocentric distance in the global mass–morphology plane. We find a significant resolved mass density–metallicity relation |$r\Sigma _\star ZR$| for galaxies of all types and masses above |$10^{9.8}\, \mathrm{M_\odot }$|⁠. Different radial distances make an important contribution to the spread of the relation. Particularly, in low- and intermediate-mass galaxies, we find that at fixed Σ_⋆ metallicity increases with radius independently of morphology. For high masses, this radial dependence is only observed in high Σ_⋆ regions of spiral galaxies. This result calls for a driver of metallicity, in addition to Σ_⋆ that promotes chemical enrichment in the outer parts of galaxies more strongly than in the inner parts. We discuss gas accretion, outflows, recycling, and radial migration as possible scenarios.

techniques: spectroscopic, galaxies: abundances, galaxies: evolution, galaxies: statistics, galaxies: stellar content

1 INTRODUCTION

Stellar metallicity saves imprints of the star formation activity and recycling efficiency in a galaxy across cosmic time and is key to our understanding of the formation and evolution of galaxies. While gas-phase metallicity reflects the current state of metal abundances and can be understood as the integral over the chemical enrichment (and dilution) history, the metallicity in stars provides evidence of its past evolution.

Numerous works have studied global galaxy parameters in the context of star formation regulation and a set of fundamental relations has been established, in particular between the total stellar mass of a galaxy and its star formation rate (SFR; star formation main sequence, SFMS; Brinchmann et al. 2004; Elbaz et al. 2007; Noeske et al. 2007; Renzini & Peng 2015), between the gas density and star formation rate (Schmidt–Kennicutt law; Schmidt 1959; Kennicutt 1998), and between the stellar mass and the metallicity, the global mass–metallicity relation (MZR) both in the interstellar medium (ISM; Tremonti et al. 2004; Kewley & Ellison 2008; Lequeux et al. 2009; Mannucci et al. 2010; Rosales-Ortega et al. 2012; Lara-López et al. 2013; Lian et al. 2015; Curti et al. 2020) and in stars (Gallazzi et al. 2005; Thomas et al. 2005, 2010; Panter et al. 2008; Kirby et al. 2013; Peng, Maiolino & Cochrane 2015; Zahid et al. 2017; Trussler et al. 2020).

With the advent of integral field spectroscopy (IFS) and the possibility of spatially resolved spectroscopic observations, the study of these fundamental relations has advanced to individual spatial regions of kpc scales within galaxies. There is growing evidence that the global relations have their local counterpart. On kpc scales, stellar surface mass density Σ_⋆ correlates with SFR (Sánchez et al. 2013; Wuyts et al. 2013; Cano-Díaz et al. 2016, 2019), metallicity (Rosales-Ortega et al. 2012; Sánchez et al. 2013; González Delgado et al. 2014b; Barrera-Ballesteros et al. 2016, 2017), and molecular gas mass (Lin et al. 2019, 2020; Ellison et al. 2020).

Spatial variations of galaxy properties have been studied for decades (e.g. Pagel & Edmunds 1981; Vila-Costas & Edmunds 1992) with increasing intensity in recent years owing to the emergence of large observational data sets of IFS surveys, such as Spectroscopic Areal Unit for Research on Optical Nebulae (SAURON; de Zeeuw et al. 2002), ATLAS^3D (Cappellari et al. 2011), Calar Alto Legacy Integral Field Area (CALIFA; Sánchez et al. 2012), Sydney-AAO Multi-object Integral field spectrograph (SAMI; Croom et al. 2012; Bryant et al. 2015), Mapping Nearby Galaxies at Apache Point Observatory (MaNGA; Bundy et al. 2015), and All-weather MUse Supernova Integral-field of Nearby Galaxies (AMUSING; Galbany et al. 2016). Research has focused particularly on gradients of metallicity in order to study the assembly history of galaxies and inflow and outflow mechanisms.

Most of the works on gas-phase abundance gradients clearly indicate negative slopes with relatively low scatter when the radius is scaled to the size of the galaxy (effective radius R_e or R₂₅; Sánchez et al. 2014; Ho et al. 2015; Sánchez-Menguiano et al. 2016). While some studies find a trend of metallicity gradient with stellar mass, in that lower mass galaxies possess flatter gradients (Belfiore et al. 2017; Poetrodjojo et al. 2018), other works suggest no dependence on stellar mass (Lian et al. 2018; Sánchez-Menguiano et al. 2018; Bresolin 2019). Other studies have suggested a trend between metallicity gradient and size (e.g. Carton et al. 2018), with Boardman et al. (2021) reporting a significant trend in gradients across the galaxy mass–size plane. The discrepancies seem to arise from different adopted methods to calibrate the abundances, different radial coverage, and, potentially, from different spectral resolutions (see e.g. Maiolino & Mannucci 2019; Sánchez 2020, for a review).

Regarding radial trends of stellar metallicity, there is a general consensus that metallicity gradients are negative in massive galaxies yet become flatter with decreasing mass (Sánchez-Blázquez et al. 2014; González Delgado et al. 2015; Goddard et al. 2017; Li et al. 2018; Lian et al. 2018; Domínguez Sánchez et al. 2019; Oyarzún et al. 2019; Lacerna et al. 2020). The correlation with mass seems to be stronger for late-type galaxies (LTGs) than for early-type galaxies (ETGs; Goddard et al. 2017).

During the last 10 yr, much observational evidence has been provided for a local relation between surface mass density and gas-phase metallicity (Rosales-Ortega et al. 2012; Sánchez et al. 2013; Barrera-Ballesteros et al. 2016), which has also been recently reproduced in cosmological simulations (e.g. Trayford & Schaye 2019). On the other hand, the spatially resolved stellar surface mass density–stellar metallicity relation (⁠|$r\Sigma _\star ZR$|⁠) has not yet been studied in great detail. It was explored in González Delgado et al. (2014b, hereafter GD14) based on a sample of 300 galaxies from the CALIFA survey, and it was shortly presented in a new collection of data from different surveys comprising 1494 galaxies in Sánchez (2020). Interestingly, Zhuang et al. (2019) find a correlation of stellar metallicity with dynamically derived surface mass density from Schwarzschild modelling of 244 CALIFA galaxies. Finally, Zibetti et al. (2020) explored the |$r\Sigma _\star ZR$| for 69 CALIFA ETGs. A more in-depth comparison to these works is presented in Section 4.

Given the connection between surface mass density and metallicity, we can partially explain a radially decreasing metallicity with an accordingly decreasing surface brightness profile (de Vaucouleurs 1959; Sersic 1968; Freeman 1970) and, thus, a decreasing surface mass density (Bakos, Trujillo & Pohlen 2008; García-Benito et al. 2019). However, it is not yet clear whether the surface mass density is sufficient in explaining radial variations of stellar metallicity or whether other additional drivers of metallicity may play a role, such as the total stellar mass, the morphological type, the environment, accretion and outflow of gas, or structural galaxy components such as galaxy bars and bulges.

In this work, we aim to investigate the distribution of stellar metallicity in nearby galaxies as a function of galactocentric radius and stellar surface mass density in the global mass–morphology plane. We make use of an unprecedented large sample of galaxies from the MaNGA survey, which allows us to explore the relations in several bins of stellar mass and morphology while keeping the statistics sufficiently large.

This paper is organized as follows. In Section 2, we give a brief outline of the MaNGA survey and our sample selection. This is followed by a presentation of our analysis in Section 3, where we introduce an update to our MaNGA firefly catalogue and derive an empirical prescription for galaxy inclination corrections. Thereafter, we present our results in Section 4, discuss our findings in Section 5, and conclude with Section 6.

Throughout the paper we assume a flat cosmology with a Hubble constant of |$H_0 = 67.8\mathrm{\, km\, s^{-1}\, Mpc^{-1}}$| and Ω_m = 0.308 (Planck Collaboration XIII 2016).

2 DATA

2.1 The MaNGA survey

This work is part of the Mapping Nearby Galaxies at Apache Point Observatory (MaNGA) survey (Bundy et al. 2015), which is a Sloan Digital Sky Survey-IV (SDSS-IV) project (Blanton et al. 2017). MaNGA uses specially designed integral field unit fibre systems (Drory et al. 2015) that feed into the Baryonic Oscillation Spectroscopic Survey (BOSS) spectrograph (Smee et al. 2013) mounted at the Sloan Foundation 2.5-m Telescope (Gunn et al. 2006).

Each observing plate employs an array of 17 hexagonal bundles of different sizes between 19 and 127 fibres to ensure a uniform spatial coverage per galaxy out to at least |$1.5\, R_\mathrm{ e}$| for the Primary+ MaNGA sample (∼2/3 of the total sample) and out to |$2.5\, R_\mathrm{ e}$| for the higher redshift Secondary sample (∼1/3 of the total sample; Law et al. 2015; Wake et al. 2017) at a median spatial resolution of |$2.54\, \mathrm{ arcsec}$| full width at half-maximum (FWHM; Law et al. 2016). Objects are observed 2–3 |$\rm h$| to reach the signal-to-noise ratio (S/N) per pixel required by the survey goals (Yan et al. 2016b). The spectra cover a wavelength range of 3622–10354 Å at a spectral resolution of |$\rm \sigma = 72\, km\, s^{-1}$| (Law et al. 2016).

Raw observations are spectrophotometrically calibrated (Yan et al. 2016a) and reduced by the Data Reduction Pipeline (drp; Law et al. 2016). High-level data products for each galaxy such as stellar kinematics, emission line properties, and spectral indices are produced by the Data Analysis Pipeline (dap; Belfiore et al. 2019; Westfall et al. 2019).

The MaNGA survey has concluded its observations in 2020 August. The final sample encompasses 10 010 unique galaxies across |${\sim} 4000\, \mathrm{deg^2}$| at a median redshift of z ∼ 0.03 and will be released by the end of 2021. The sample was selected to have a flat stellar mass distribution in logarithmic space in the range |$5\times 10^8 \le M_{\star } \le 3\times 10^{11}\, \mathrm{ M}_{\odot }$|⁠. Additionally, a colour-enhanced supplement sample was added to upweight underrepresented galaxies in the colour–magnitude diagram. The full main MaNGA sample consists of the Primary sample (⁠|${\sim} 50{{\ \rm per\ cent}}$|⁠) and its Colour-Enhanced supplement (⁠|${\sim} 17{{\ \rm per\ cent}}$|⁠; together called Primary+), with radial coverage out to |$1.5\, R_\mathrm{ e}$|⁠, and the Secondary sample (⁠|${\sim} 33{{\ \rm per\ cent}}$|⁠), with radial coverage out to |$2.5\, R_\mathrm{ e}$| (for more details, see Wake et al. 2017).

2.2 Sample selection

The parent sample of this work is drawn from the final MaNGA Product Launch 11 (MPL-11), which will be released in the upcoming SDSS Data Release 17 (DR17). The data comprise 11 273 data cubes, 10 145 of which represent high-quality galaxy observations corresponding to 10 010 unique galaxies and 135 repeat observations. From the 10 010 galaxies, we selected all those that are part of the main MaNGA sample¹ and that have been analysed by the drp and dap, which yields a sample of 9714 galaxies.

From this sample, we removed

• data cubes that were classified as critical by the drp or by the dap;

• galaxies from the galaxy pair ancillary programmes;

• galaxies with inclinations i > 80° (see Section 3.2);

• data cubes with less than 30 Voronoi bins for a given galaxy (to assure decent spatial sampling, see Section 3.1).

Our analysis further requires a morphological classification. For that purpose we cross-matched our sample with that of the MaNGA Deep Learning Morphological Value-Added-Catalogue-Data Release 15 (MDLM-VAC-DR15; Fischer, Domínguez Sánchez & Bernardi 2019) that uses Deep Learning algorithms based on Convolutional Neural Networks described in more detail in Domínguez Sánchez et al. (2018). In this work, we use an updated version of the MDLM-VAC that contains 10 293 unique entries from the upcoming DR17, kindly provided by H. Domínguez Sánchez and collaborators (MDLM-VAC-DR17; Domínguez Sánchez et al., in preparation). Morphological T-types are obtained from models that have been trained and tested with the morphological galaxy catalogue from Nair & Abraham (2010) on Red–Green–Blue (RGB) SDSS Data Release 7 (DR7) cutout images. Values of the T-types range from −4 to 9. Two additional classifications from the catalogue are used in this paper: P_S0 for a finer separation between ellipticals (Es) and lenticulars (S0), and P_LTG to help separating ETGs from LTGs in addition to the T-types. The value given in the catalogue is the average of the output of 10 different models trained with k-folding. Finally, we also include the visual classification flag VC. We use the following definition to separate between Es (T ≤ 0, P_LTG < 0.5, P_S0 < 0.5, and VC = 1), S0 (T ≤ 0, P_LTG < 0.5, P_S0 ≥ 0.5, and VC = 2), and LTGs (T > 0, P_LTG ≥ 0.5, and VC = 3). The cross-match of our selected sample with the objects in the MDLM-VAC that have a morphological classification with the aforementioned criteria yields our final working sample of 7439 galaxies (2215 Es, 808 S0s, and 4416 LTGs).

Fig. 1 presents the sample in the near-ultraviolet (NUV) − r colour–mass and redshift–mass diagram in comparison with the main DR17 sample. We use total galaxy masses from the NASA-Sloan Atlas (NSA)² that were obtained from K-correction fits to elliptical Petrosian photometric fluxes (Blanton et al. 2011; Blanton & Roweis 2017), assuming a Chabrier (2003) initial mass function (IMF) and Bruzual & Charlot (2003) single stellar population models (SSPs). These fits are held to contain nearly all flux of a galaxy and do not depend on aperture correction. We use the total galaxy masses in this work for the global relations shown in Section 4.1 and to separate the sample in stellar mass bins. We scaled the masses to the Hubble parameter |$H_0 = 67.8\mathrm{\, km\, s^{-1}\, Mpc^{-1}}$| used in this work.

Figure 1.

Top: colour–mass diagram derived from SDSS photometric masses, NUV- and r-band filters, taken from the NSA catalogue. Grey contours and histograms show the complete MaNGA DR17 main sample. Scatter points show individual galaxies from the selected working sample separated by morphology into ellipticals (Es, red), lenticulars (S0, orange), and late-type galaxies (LTG, blue). Histograms for E, S0, and LTG are stacked. Bottom: redshift–mass distribution of the sample. Redshifts are taken from the spectroscopical estimate ‘stellar_z’ in the dap.

There is no obvious bias in the sample and our selection seems representative of the complete MaNGA sample. The distribution is close to flat over a large range of mass and colour. Galaxies classified in the MDLM-VAC as LTGs occupy predominantly the blue cloud, whereas S0s and Es are located on the red sequence and dominate the high-mass end of the distribution. The gap seen in the redshift distribution is a natural consequence of the higher redshift selection criterion of the Secondary sample as compared to the Primary+ sample.

3 ANALYSIS

The results presented in this work are based on measurements of spatially resolved stellar population properties, such as stellar metallicities ([Z/H]), stellar ages, masses (M_⋆), and surface mass densities (Σ_⋆). The derivation of these parameters is a multistep process that is outlined in the following.

3.1 MaNGA firefly VAC

The extraction of the stellar population properties is done by employing the full spectral fitting technique using our own firefly³ code (Wilkinson et al. 2017). firefly is a χ² minimization code that fits combinations of SSPs to spectroscopic data.

The firefly analysis of MaNGA data builds upon the dap, which already includes all necessary preparatory steps. In particular, the dap workflow include (1) adaptive Voronoi spatial binning (Cappellari & Copin 2003) to a minimum target |$\rm S/N \sim 10$| in the stellar continuum, (2) measurements of stellar kinematics, and (3) emission line fluxes employing the full spectral fitting code ppxf (Cappellari & Emsellem 2004; Cappellari 2017). The dap uses a hierarchically clustered (HC) set of template spectra from the MILES stellar library (Sánchez-Blázquez et al. 2006) for the extraction of stellar kinematics, while the stellar continuum modelling in the emission line module uses a subset of MaStar SSPs (Maraston et al. 2020) based on the MaStar stellar library (Yan et al. 2019).

From each dap-processed data cube, firefly reads the Voronoi-binned spectra, subtracts emission line fluxes, and shifts the spectra to rest-frame wavelengths using the redshift and stellar velocity information from the dap. Additionally, the stellar velocity dispersion and the instrumental resolution are used to match the line broadening of the SSPs to the data prior to the fitting process. For the results presented in this work, we used the stellar population models of Maraston & Strömbäck (2011) based on the MILES stellar library (m11-MILES) and a Kroupa (2001) IMF. Updates and comparison using our new MaStar models (Maraston et al. 2020) are subject for future papers.

firefly fits arbitrarily weighted linear combinations of SSPs to the data following an iterative process controlled by the Bayesian information criterion. No additive or multiplicative polynomials are employed and star formation histories (SFHs) are not regularized in order to allow for a large fitting freedom and sufficient exploration of the parameter space. firefly has been shown to perform well down to |$\rm S/N \sim 5$| (see Wilkinson et al. 2017, for more detail).

The results of the complete firefly run on all DR17 data cubes are stored in the MaNGA firefly Value-Added-Catalogue (VAC; Neumann et al., in preparation).⁴ This constitutes an update of the earlier versions of the VAC (Goddard et al. 2017; Parikh et al. 2018) and provides spatially resolved stellar population properties, including stellar mass and its partition into stellar remnants, age, metallicity, stellar surface mass density, dust attenuation, and the full SFHs for 10 735 data cubes. The MaNGA firefly VAC complements the dap by providing further high-level data products in addition to the stellar kinematics, emission line properties, and spectral indices given by the dap.

Throughout the analysis in this paper, we use light-weighted and linearly averaged stellar population parameters. In most cases, however, we employ the median, which is identical in linear and logarithmic space. In addition to our galaxy sample selection, we apply a |$\rm S/N\gt 8$| cut to all Voronoi bins. Furthermore, we correct all stellar ages to the median redshift of the sample z = 0.0376.

3.2 Inclination correction

The stellar surface mass density Σ_⋆ is defined as the stellar mass divided by the surface area. The stellar masses are calculated by firefly for each spatial Voronoi bin and Σ_⋆ is given in the MaNGA firefly VAC, but the area of the corresponding bin still needs to be corrected for projection effects. For a circular thin disc the deprojected area is given by A = A_obs (cos i)⁻¹, where A_obs is the observed area and i is the inclination of the disc. This yields |$\Sigma _\star = \Sigma _{\star ,\, \mathrm{obs}} \, \cos i$| for the surface mass density. However, real galaxies have an intrinsic thickness q₀, the exact value of which is often very uncertain and varies with morphology (Heidmann, Heidmann & de Vaucouleurs 1972; Bottinelli et al. 1983; Guthrie 1992), wavelength (Mitronova et al. 2004), and kinematics (Weijmans et al. 2014). The inclination of a thick oblate spheroid is obtained as |$\cos i = \sqrt{(q^2-q_0^2)/(1-q_0^2)}$|⁠, where q is the observed axial ratio of the projected spheroid (Hubble 1926). Thus, we have

$$\begin{eqnarray*} \Sigma _\star = \Sigma _{\star ,\, \mathrm{obs}} \, \sqrt{\left(q^2-q_0^2 \right)/ \left(1-q_0^2 \right)}. \end{eqnarray*}$$

(1)

The apparent axial ratios of the galaxies in our sample are obtained from the elliptical Petrosian analysis in the enhanced NSA catalogue (Wake et al. 2017).

Under the simplification that the intrinsic axial ratios only vary with morphology and that galaxies seen face-on are perfectly circular (but see e.g. Ryden 2004, who find that discs are slightly elliptical), we derive the intrinsic thickness of each galaxy following a similar procedure as described in Heidmann et al. (1972), Bottinelli et al. (1983), and Guthrie (1992).

In Fig. 2, we plot the observed axial ratio of each galaxy in the sample against the morphological T-type. In this approach, the assumption is that at a given T-type the galaxies with the lowest axial ratio are seen edge-on and, thus, for these galaxies q₀ = q. We fit a linear function in log-linear space to the lower 2.5th percentile of the sample to allow for scatter due to measurement errors, as seen in the figure. This yields the relation

$$\begin{eqnarray*} {-}\log q_0 = 0.316 + 0.049\, T. \end{eqnarray*}$$

(2)

Figure 2.

Axial ratios q = b/a as a function of morphological T-type. Black points are individual galaxies. The dark grey line shows the running median and the grey shaded areas from dark to light mark the inner |$68{{\ \rm per\ cent}}$|⁠, |$95{{\ \rm per\ cent}}$|⁠, and |$99{{\ \rm per\ cent}}$| around the median, respectively. The blue dashed line is a linear fit to the 2.5th percentile.

We tested whether the derivation of q₀ is dependent on galaxy mass by selecting galaxies of different mass bins and reevaluating Fig. 2 and found no significant change. This equation is directly empirically derived from the same sample of MaNGA galaxies for that we aim to derive the inclinations in this work, but we would like to emphasize that it can be used for any kind of galaxy given the assumptions stated above.

In comparison with our analysis, we also show in Fig. 2 the relations derived by Heidmann et al. (1972, using the measurements from Holmberg 1958), Bottinelli et al. (1983), Paturel et al. (1997), and Bizyaev et al. (2014). The slope of our equation (0.049) is very similar to those in the literature (0.047, 0.045, 0.053, and 0.0497, respectively). However, the derivation of the intercept is clearly different, with our measurements yielding the largest axis ratios of edge-on systems. Part of the disagreement stems from the somewhat arbitrary choice of allowing for a certain amount of scatter above the fitted upper limit. The more recent results from the mentioned literature, Paturel et al. (1997) and Bizyaev et al. (2014), for example, agree well with our statistics, would we be more stringent with the scatter. In addition to that, the measurement of the axis ratios of galaxies is sensitive to the applied technique. For instance, the authors cited above use ellipse fitting at a surface brightness level of 25 and 26.5 |$\rm mag\, arcsec^{-2}$| or at a certain S/N level, while the axis ratios from the NSA catalogue employed here are derived at the 90 per cent Petrosian light radius. Furthermore, some authors who find very low q₀ ≤ 0.2 for LTGs correct axial ratios for the distortions of the bulge component (cf. Giovanelli et al. 1994) what naturally leads to smaller ratios for smaller T-types.

Comparing our results with the kinematic analysis of ETGs in Weijmans et al. (2014), we find a very good agreement for our galaxies with T ≤ 0 and 0.75 > q₀ > 0.48 and their results for slow rotators: q₀ = 0.63 ± 0.09. Fast rotators, on the other hand, have been found to be much flatter with q₀ = 0.25 ± 0.14, a range of axis ratios that is only populated by LTGs in Fig. 2.

Using equations (1) and (2), we are now able to determine the inclination and intrinsic surface mass density. The distribution of inclinations is shown in Fig. 3. To further test our results we compared a subset of ∼450 galaxies to the kinematically derived inclinations for the same galaxies in Yang et al. (2021) and found a good agreement.

Figure 3.

Distribution of galaxy inclinations derived with equation (2) from Fig. 2.

In addition to conducting an inclination correction to the surface mass density, we decided to exclude galaxies with i > 80° from the further analysis, because edge-on galaxies are prone to large uncertainties in the derivation of integrated stellar population measurements. We tested how different inclination cuts affected our analysis and found that the main results of this work remain unchanged, but including galaxies with i > 80° introduced increased scatter into the derived relations.

4 RESULTS

In this section, we present the results of our analysis. We start by reproducing the well-known global MZR, but separated by morphology, and discuss it in light of previous studies. Afterwards, we zoom-in on spatially resolved relations starting with radial variations of stellar metallicity before presenting the complex interplay between metallicity, surface mass density, galactocentric distance, total mass, and morphology from different angles.

4.1 Global stellar mass–age and mass–metallicity relation

The global mass–metallicity relation (MZR) is a relation between the total stellar mass of a galaxy and its metallicity. For this relation, measurements of the central metallicity (e.g. Trussler et al. 2020), the metallicity at the effective radius R_e, the galaxy-wide averaged metallicity (GD14), or the mean metallicity within R_e (Sánchez et al. 2018) have been used in the literature. In all cases a clear correlation has been found, indicative of an underlying local relation.

In Fig. 4, we show the global mass–age and mass–metallicity relation of ETGs and LTGs in our sample. Both age and metallicity are averaged within the central 3 arcsec diameter for a better comparison with previous SDSS-based literature. These values are provided in our MaNGA firefly VAC. We do not find a significant difference regarding the central stellar population properties between elliptical and lenticular galaxies and decided to combine them here as ETGs. They will be shown separately in the rest of the work.

Figure 4.

Global mass–age relation (top) and mass–metallicity relation (bottom) separated by morphology. Age and metallicity are averaged within a central area of 3 arcsec diameter. The scatter points represent individual galaxies and the line plot shows the median age and metallicity across mass bins of |$0.2$| dex width. The error bars display the median absolute deviation.

The global mass–age relation for LTGs is mildly positive for |$M_\star \lt 10^{10}\, \mathrm{M_\odot }$| with a stronger increase towards higher masses. The ages of ETGs are systematically older and the distribution across different masses is also relatively flat with a slightly positive age–mass correlation. These trends are in good agreement with results from Trussler et al. (2020) who studied the global chemical properties of ∼80 000 local galaxies in SDSS, separated by SFR into passive, star-forming and green valley galaxies. In this comparison we assume that the morphological classification used in our work is a good first-degree approximation of star-forming and passive galaxies in Trussler et al. (2020). It is important to point out that while the trends are similar, the absolute values differ owing to the fact that we show light-weighted quantities, while they show mass-weighted values, the latter of which are usually older and more metal poor.

We observe a clearly positive mass–metallicity relation of both ETGs and LTGs with very similar slopes, but ETGs are systematically offset towards higher metallicities. This is again in agreement with Trussler et al. (2020) and has been explained by the authors as due to a combination of an extended starvation phase and outflows (see also Peng et al. 2015). In short, the star-forming progenitor of a passive galaxies stops accreting new gas during the quenching phase and consequently uses up the remaining gas in the ISM efficiently and without dilution to form new metal-rich stars, which leads to a sharp increase in stellar metallicity while stellar mass increases only slowly.

In summary, we find a mildly positive global mass–age relation and a positive mass–metallicity relation, with the central parts of ETGs being systematically older and more metal rich than LTGs at fixed total mass. The positive MZR from our IFS data reproduces well earlier findings from both single-fibre observations (e.g. Thomas et al. 2010; Zahid et al. 2017; Trussler et al. 2020), semi-analytic models (Henriques et al. 2020; Yates et al. 2021), and cosmological simulations (e.g. Schaye et al. 2015; Ma et al. 2016; De Rossi et al. 2017).

4.2 Radial trends of stellar populations

After testing the global galaxy-to-galaxy variations of stellar populations with stellar mass and morphology, we are now in the position to further exploit the potential of the IFS data and study stellar population parameters locally as a function of galactocentric distance. Henceforth, we separate our galaxy sample in three morphology bins (elliptical, lenticular, and spiral) times five total stellar mass bins (log M_⋆/M_⊙ ∈ [7.3, 9.8, 10.3, 10.8, 11.3, 12.3]). We further average the spatially resolved parameters (from the Voronoi bins) of each galaxy individually in radial bins of 0.05 × R_e size. This is done to prevent galaxies with many Voronoi bins from being upweighted and, thus, in order to assure that each galaxy contributes equally to the subsequent analysis. In the following, when we talk about individual spatial bins, we are not referring to the Voronoi bins but to the radial bins. Examples of 2D stellar population maps produced by firefly from nine galaxies across the mass–morphology grid can be found in Fig. A1.

4.2.1 Age and metallicity profiles

In Fig. 5, we present radial age and metallicity profiles of our sample. The coloured density plots in the background show the distribution of individual spatial bins, while the black line shows the median in radial bins of |$0.1\, R_\mathrm{e}$|⁠. Additionally, we show a linear fit to the profile within |$1\, R_\mathrm{e}$|⁠, which is plotted as solid orange line for R ≤ R_e and as dashed extrapolated line for larger radii.

Figure 5.

Radial metallicity (top) and age (bottom) profiles. Columns show different bins of total stellar mass M_⋆. Rows show different morphologies. Each panel presents a density plot of individual spatial bins from galaxies with the corresponding mass and morphology. The data are smoothed with a Gaussian kernel density estimator. In black, we show median [Z/H] and age in bins of |$0.1\, R_\mathrm{e}$| width. The error bar is the standard error on the median, calculated as |$\sigma _\mathrm{err}=(\pi /2)\, (\sigma /\sqrt{N})$|⁠, where σ is the standard deviation and N the sample size per bin. A linear regression to the data within |$1\, R_\mathrm{e}$| is plotted as solid orange line with an extrapolation towards larger radii as dashed line. The corresponding slope is annotated in each panel.

We observe negative inner metallicity gradients for all galaxies with log M_⋆/M_⊙ > 9.8, which flatten towards larger radii, most clearly seen for Es and intermediate-mass S0s. For the lowest mass galaxies, metallicity profiles of Es and S0s are very shallow across the whole radial range and are clearly positive (⁠|$\nabla [Z/H]_{1\, R_\mathrm{e}} = 0.057 \pm 0.004\, \mathrm{dex}/R_\mathrm{ e}$|⁠) for spiral galaxies. Furthermore, it is interesting to note that while the metallicity changes dramatically in the centre across the mass–morphology plane – giving rise to the global MZR as presented in the previous section – the outermost parts at |$R=2.5\, R_\mathrm{e}$| do not change much. This could be an indication of a relatively uniform stellar metallicity in the outer halo of galaxies. However, we do not sample these regions well enough to draw firm conclusions.

The age profiles show relatively flat, mildly positive gradients for elliptical and lenticular galaxies for all masses (⁠|$|\nabla \mathrm{Age}_{1\, R_\mathrm{ e}}|\lt 1\, \mathrm{Gyr}/R_\mathrm{ e}$|⁠). They are negative but shallow for spiral galaxies with a central steepening at high masses (log M_⋆/M_⊙ > 10.8), indicative of an old bulge component. Nevertheless, it is apparent in the density plot that not all LTGs in the high-mass bin have this central upturn. In fact, a substantial fraction seems to have a rather flat age distribution pointing towards the absence of a dominating classical bulge. Lenticular galaxies show an interesting transition in the age distribution from mostly older ages for log M_⋆/M_⊙ < 11.3 to younger ages for log M_⋆/M_⊙ > 10.8, with a bimodality in between which is typically observed between active and passive galaxy populations. This could be related to the differences seen in the stellar population gradients between low- and high-mass S0s in Domínguez Sánchez et al. (2020), but needs yet to be investigated in more detail.

Increasing the sample size of previous works by more than an order of magnitude, we find that the median age and metallicity gradients presented here in light-weighted quantities and their trends with mass and morphology are qualitatively in good agreement with previous analysis of MaNGA data with the firefly full spectral fitting code (Goddard et al. 2017) and other literature (e.g. Kuntschner et al. 2010; Sánchez-Blázquez et al. 2014; González Delgado et al. 2015; Li et al. 2018; Lian et al. 2018; Oyarzún et al. 2019; Zhuang et al. 2019; Lacerna et al. 2020; Sánchez 2020). However, the metallicity gradients for massive Es (⁠|$\nabla [Z/H]_{1\, R_\mathrm{e}} = -0.116 \pm 0.002\, \mathrm{dex}/R_\mathrm{ e}$|⁠) and S0 (⁠|$\nabla [Z/H]_{1\, R_\mathrm{e}} = -0.073 \pm 0.012\, \mathrm{dex}/R_\mathrm{ e}$|⁠) are flatter than some results from the literature indicate (e.g. Kuntschner et al. 2010; Martín-Navarro et al. 2018; Domínguez Sánchez et al. 2019; Zibetti et al. 2020), yet they are consistent with other studies (e.g. Mehlert et al. 2003; González Delgado et al. 2015; Goddard et al. 2017; Li et al. 2018; Lacerna et al. 2020). Some of the discrepancy can be attributed to differences in stellar population analysis techniques (full spectral fitting, line strength measurement), averaging (light-weighted, mass-weighted, linear, logarithmic), sample selection, and data quality (S/N, sampling); see e.g. Goddard et al. (2017) and Zibetti et al. (2020) for further discussion. The detailed variations seen in the density distributions of single spatial bins across the mass–morphology plane leave room for further research in upcoming papers.

4.2.2 Surface mass density profiles

In this work, we aim to relate these stellar population gradients with stellar mass surface density for which we present radial profiles in Fig. 6. For simplicity we only show the medians for each mass–morphology bin. From this figure, it is apparent that Σ_⋆ is monotonically decreasing with radius for all galaxies. Higher mass galaxies have higher Σ_⋆ at any given radius (except for the highest mass ETGs) and a more pronounced central upturn. Furthermore, spiral galaxies have lower Σ_⋆ at fixed mass and radius than S0 or Es. The trends are expected to change if the radius is not normalized by the effective light radius R_e due to the correlation between R_e and M_⋆, as well as its dependence on morphology (e.g. van de Sande et al. 2019; Boardman et al. 2021). These results are qualitatively consistent with previous studies (e.g. González Delgado et al. 2014a; García-Benito et al. 2017; Zibetti et al. 2020).

Figure 6.

Radial profiles of stellar surface mass density in different bins of total stellar mass M_⋆ (grey scale) and morphology (panels from top to bottom). Right-hand panels show the equivalent profiles from the left-hand side but normalized to the central surface mass density Σ_{⋆, 0}.

The exclusively negative Σ_⋆ gradients and the mostly negative metallicity gradients suggest that a local correlation should be found between Σ_⋆ and [Z/H] as local spatially resolved counterpart to the MZR, the |$r\Sigma _\star ZR$|⁠. At the same time, however, the differences seen in the gradients already point towards additional radially dependent drivers of stellar metallicity.

4.3 Local stellar surface mass density–metallicity relation

The local and spatially resolved stellar surface mass density–stellar metallicity relation, |$r\Sigma _\star ZR$|⁠,⁵ is shown in the left-hand panel of Fig. 7. This is a density plot based on >2.6 million Voronoi bins, radially rebinned as described above, from the complete working sample including galaxies from all morphologies, masses, and from all radial positions out to |$2.5\, R_\mathrm{e}$|⁠.

Figure 7.

Spatially resolved stellar surface mass density–metallicity relation. Left: density plot of all spatial bins from the complete sample out to maximal |$2.5\, R_\mathrm{e}$|⁠. The data are smoothed with a Gaussian kernel density estimator. Contours enclose |$20{{\ \rm per\ cent}}$|⁠, |$40{{\ \rm per\ cent}}$|⁠, |$60{{\ \rm per\ cent}}$|⁠, and |$80{{\ \rm per\ cent}}$| of the sample. In white, we show median [Z/H] in log Σ_⋆ bins of |$0.1$| dex width in the range 7.3 < log Σ_⋆ < 9.4 and, in green, a linear regression to the data within the same range. The standard error on the median, calculated as |$\sigma _\mathrm{err}=(\pi /2)\, (\sigma /\sqrt{N})$|⁠, where σ is the standard deviation and N the sample size per bin, is smaller than the marker size and thus not visible. The Spearman’s rank correlation coefficient ρ is shown in the bottom right. Right: same as left, but coloured according to the median age in each hexagonal bin without smoothing. Only bins with more than 20 data points are shown.

The figure shows a positive log Σ_⋆–[Z/H] relation with a monotonic and close to linear increase of the binned median values. Because of a relatively broad distribution, we find a mild Spearman’s rank correlation coefficient of ρ = 0.64 with a p-value of p ≪ 0.001. A linear regression yields the general equation:

$$\begin{eqnarray*} [Z/H] = (0.221 \pm 0.001)\, \log \Sigma _\star - (1.803 \pm 0.006). \end{eqnarray*}$$

(3)

In the right-hand panel of Fig. 7, we display the same plot of the |$r\Sigma _\star ZR$| but colour coded according to the median age. A clear positive trend of stellar age with both stellar metallicity and surface mass density is visible. However, the oldest and youngest population are not found at the highest and lowest metallicities but rather at the densest and less densest regions, respectively. Yet, at fixed mass density older populations tend to have higher metallicity.

From what we have seen, it seems clear that surface mass density drives stellar population properties locally, but we would also like to understand what causes the spread in the distribution and whether additional trends can be identified. A first step in this direction is made in Fig. 8, where we present once again the |$r\Sigma _\star ZR$|⁠, this time coloured by the median galactocentric distance. At fixed metallicity, it is apparent that surface mass density is a decreasing function of radius corroborating the relations illustrated in Fig. 6. Interestingly, at fixed Σ_⋆, metallicity increases with radius. In fact, following same coloured regions, stellar populations at the same median radius display a relatively tight |$r\Sigma _\star ZR$|⁠. This relation is shifted towards higher metallicities and lower surface mass densities at larger radii.

Figure 8.

Same as Fig. 7, but coloured according to the median radius in each hexagonal bin without smoothing. Only bins with more than 20 data points are shown. Contours show smoothed |$80{{\ \rm per\ cent}}$| density distributions from the inner part of galaxies, with |$R \le 0.25\, R_\mathrm{e}$|⁠, in different bins of total stellar mass M_⋆.

Furthermore, we explore how the total mass of the host galaxy affects the relation at fixed radius. In Fig. 8, we superimpose contours for different mass bins and we consider only the central regions (⁠|$R\le 0.25\, R_\mathrm{e}$|⁠). Each contour encloses 80 per cent of the subsample. The relation between [Z/H] and Σ_⋆ is reproduced at each mass bin (although it is not clear whether for the lowest mass bin we observe a very steep correlation or no correlation at all). As total mass increases, the distribution moves along a slightly curved |$r\Sigma _\star ZR$| with flattening slope towards higher masses. This behaviour is also representative for the trends at larger (but fixed) radii, which we not show here to maintain the clarity of the figure.

The dependence on total galaxy stellar mass and morphology looks different, if we do not separate the relation by radius but instead plot galaxy-wide |$r\Sigma _\star ZR$| separated by mass and morphology as shown in Fig. 9. In contrast to the mass contours in Fig. 8, this plot shows all spatial bins across each field of view (FoV); a mass–morphology separated version of Fig. 7.

Figure 9.

Spatially resolved stellar surface mass density–metallicity relation. Same as left-hand panel of Fig. 7, but here separated by total stellar mass and morphology of the respective galaxy. Sky blue dots show the median [Z/H] in |$\log \, \Sigma _\star$| bins of 0.1 dex.

A positive rΣ_⋆ZR is found for all masses and morphologies except for the lowest mass bin. However, the correlation becomes weaker for lower galaxy masses and lower local mass densities. In fact, we observe an interesting bending in the relation between Σ_⋆ ∼ 10⁸ and |$10^9\, \mathrm{M_\odot \, kpc^{-2}}$| at lower galaxy masses and almost independent of morphology. Similar upturns at high surface mass densities are also apparent in Sánchez (2020, fig. 9).

Thus, from Fig. 9 we conclude that the shape of the rΣ_⋆ZR depends on total galaxy mass and morphology if radius is not fixed. This is qualitatively in good agreement with the results presented in GD14 and Sánchez (2020). We speculate that the disappearance of a mass density dependence at low total masses and low local mass densities might be an effect of increasing and outweighing importance of an additional radial-dependent driver of metallicity.

Altogether, the spatially resolved stellar population properties presented in Figs 7–9 paint a picture in which old and metal-rich populations are predominantly found in the centres of massive galaxies, the most metal-poor populations are found in the centres of the least massive galaxies, while the youngest populations are found in the outer parts of low-mass galaxies. The |$r\Sigma _\star ZR$| is reproduced at all masses except for the lowest mass bin. There is some indication that the relation is shifted towards higher metallicities at larger radii, but this cannot be clearly disentangled from an overlaying mass effect in Fig. 8.

Very similar trends of age, M_⋆, and R in the Σ_⋆–[Z/H] plane are found in GD14 in an analysis of 300 CALIFA galaxies. However, the scatter in the |$r\Sigma _\star ZR$| is explained in GD14 by the impact of the global metallicity driver M_⋆. While we come to the same conclusion that high-metallicity samples come from massive galaxies, we point out that radius is an important parameter in that relation and that the scatter in the |$r\Sigma _\star ZR$| is significantly reduced if the relation is shown at fixed radius. The flattening of the relation in the inner regions of high-mass galaxies is in very good agreement with GD14 and is interpreted by these authors as due to a dominant role of M_⋆ as driver of metallicity in spheroids as opposed to a Σ_⋆-dominated regime in discs. Contrary to that interpretation, Zibetti et al. (2020) show that metallicity in massive ETGs is primarily driven by surface mass density and only secondarily modulated by galaxy mass. This is attributed by the authors to steeper metallicity gradients in their work. When we compare the |$r\Sigma _\star ZR$| in Fig. 7 to the median relation for ETGs in Zibetti et al. (2020), we find a remarkably good agreement in the overlapping range of 8 < log Σ_⋆ < 9.4 with metallicity increasing from slightly subsolar to [Z/H] ∼ 0.22. However, if we consider only ETGs with |$\log M_\star /\mathrm{ M}_\odot \gt 10.3$| in Fig. 9, our values are shifted towards slightly higher metallicities.

In Zhuang et al. (2019), the authors report a tight local relation between dynamical surface mass density and stellar metallicity within |$1\, R_\mathrm{e}$| from measurements of 244 CALIFA galaxies across all morphological types. This dynamical relation is consistent with our |$r\Sigma _\star ZR$|⁠, in particular for the central regions, and strengthen further our results.

4.4 Metallicity trends at fixed surface mass density

In Fig. 10, we present the local relations between Σ_⋆, [Z/H], and R/R_e from a different perspective. In this plot, we focus on how the metallicity changes as a function of radius at fixed surface mass density in the total stellar mass–morphology plane. We thereby test in more detail the hypothesis drawn from Fig. 8 that metallicity increases with radius at fixed mass density and fixed total mass. We plot again single spatial bins from the whole population of galaxies divided by mass and morphology. The colour of each hexagonal area shows the median [Z/H] of all contributing spatial bins in the corresponding Σ_⋆–R/R_e intervals.

Figure 10.

Stellar surface mass density versus radius coloured by median stellar metallicity for individual spatial bins. The grid shows the global mass–morphology plane of the host galaxies. The number of galaxies that contributed to each panel is shown in the lower left-hand corner. We use individual colour bar limits in each panel to highlight subtle trends. The minimum and maximum of each colour bar is shown in the lower right-hand corner of the corresponding panel. The data are smoothed using the locally weighted linear regression (loess) method as described in the main body of the text.

We employ locally weighted linear regression (loess; Cleveland & Devlin 1988) in the python implementation from Cappellari et al. (2013) to smooth the data in two dimensions. This method has been successfully tested in several works to study and recover mean trends of stellar population properties in galaxies (e.g. McDermid et al. 2015; Rosado-Belza et al. 2020; Boardman et al. 2021).

First of all, at almost any radial location in any part of the grid, we clearly see increasing metallicity with increasing surface mass density, i.e. an rΣ_⋆ZR is found at any galactocentric distance. Additionally and importantly, if we look at these plots at fixed Σ_⋆ instead of at fixed radius, we observe that for most of the mass densities, galaxy masses, and morphologies [Z/H] is increasing with radius. The radial dependence is clearest for spiral galaxies and in low-mass ETGs (log M/M_⊙ < 10.8) except for the lowest surface mass densities. It disappears in higher mass ETGs.

We explore the strength and significance of this radial dependence further in Fig. 11. For each panel across the mass–morphology grid, we compute the Spearman’s rank correlation coefficient between [Z/H] and R/R_e in fixed |$0.2\, \mathrm{dex}$| bins of log Σ_⋆. We show only statistically significant correlations by limiting the corresponding p-value to p < 0.05. Positive but weak correlations are found for K = 2, 3, 6, 7, 8, 11, 12, and 13 and at high mass densities in K = 14 and 15. This corresponds to all galaxies below log M_⋆/M_⊙ < 10.8 except for the lowest mass Es, as well as to high-density regions in massive spirals. The dependence is slightly negative in high-mass ETGs and at low-density regions in high-mass spirals. These statistics corroborate the average trends observed in Fig. 10.

Figure 11.

Spearman’s rank correlation coefficients ρ testing a correlation between [Z/H] and R/R_e in fixed bins of Σ_⋆. On the x-axis, K is a simple numbering of the 15 panels in the mass–morphology grid with K = 1–5 for ellipticals (Es) with increasing mass, K = 6–10 for lenticulars (S0), and K = 11–15 for spirals. ρ is only shown if the p-value p < 0.05 and if the number of data points is larger than 20.

In addition to local trends, we can also confirm in Fig. 10 that metallicity is globally driven by galaxy mass and morphology, in that [Z/H] is larger for more massive galaxies (from left to right) and for earlier types (from bottom to top) even at fixed Σ_⋆ and galactocentric distance. However, the cross-comparison is admittedly not an easy exercise to do in this representation given the individual colour bar limits and it is better represented in the previous sections.

From the representation of trends of stellar metallicity in Fig. 10, it is clear that negative metallicity gradients in local galaxies are largely imposed by radially monotonically decreasing surface mass density and, in fact, metallicity is not decreasing ‘steep enough’ from what mass densities would drive. At fixed surface mass density metallicity is increasing with galactocentric distance. This is true for most low- to intermediate-mass galaxies and at high Σ_⋆ in massive spirals. This result calls for additional driver of stellar metallicity that enhances metallicity more strongly at larger radii or/and dilutes it at smaller radii.

5 DISCUSSION

The observed metallicity in stars is a result of the chemical evolution and SFH. To increase the metallicity, the ejected gas from a population of stars needs to be recycled into a new stellar generation. The more generations of stars are formed the higher the possibility to build-up metallicity, but more importantly, the pre-enriched gas needs to be efficiently recycled. The relation between mass and metallicity is usually understood as a consequence of the ability to retain enriched gas ejecta. The higher the mass, the higher is the potential well that keeps gas back and efficiently recycles it to new stars (e.g. Thomas, Greggio & Bender 1999, and references therein).

Reality is, however, more complex than this closed-box scenario and cosmological simulations have shown the importance of inflows, outflows, mergers, and radial migration in shaping the radial distribution of stellar populations in galaxies (e.g. Sellwood & Binney 2002; Roškar et al. 2008; Minchev & Famaey 2010; El-Badry et al. 2016; Grand et al. 2019; Collacchioni et al. 2020).

Inverted or positive metallicity gradients at higher redshifts have been explained by a dilution of the ISM in the centre due to pristine metal-poor gas inflow in the central regions (Cresci et al. 2010) and, additionally, outflows of metal-rich gas from the centre (Troncoso et al. 2014). In the galactic fountain scenario this metal-rich gas is reaccreted on to the galaxy out to large radii (e.g. Oppenheimer et al. 2010). Higher mass galaxies are able to retain outflows earlier in cosmic history and thereby enable a rapid metal enrichment (Tremonti et al. 2004; Muratov et al. 2017). Furthermore, central metallicity dilution due to inflows can also be merger induced as shown from cosmological galaxy formation simulations in Bustamante et al. (2018), thus, environment is likely to impact the chemical evolution.

A number of cosmological simulations have shown that these and similar processes can shape the observed gas and stellar metallicity gradients in local galaxies. From the zoom-in Auriga simulations, Grand et al. (2019) report that most of the material in redshift zero stars have been ejected and reaccreted in form of galactic fountains at least once (cf. Brook et al. 2014; Übler et al. 2014; Christensen et al. 2016). By comparing different metal loading factors, they show that galactic fountains play an essential role in efficiently redistributing metal-enriched gas and thereby in flattening of the stellar metallicity gradient. However, the fractional contribution of different accretion channels is heavily debated in the literature. For example, Mitchell, Schaye & Bower (2020) find that first infall accretion clearly dominate recycled inflow in the Evolution and Assembly of GaLaxies and their Environments (EAGLE) simulations, while Anglés-Alcázar et al. (2017) report fractions for the Feedback In Realistic Environments (FIRE) simulations in between that of Auriga and EAGLE.

At later times in cosmic history metal-rich gas accretion, e.g. from intergalactic transfer (Anglés-Alcázar et al. 2017), on to the outer parts of galaxies, triggering or maintaining star formation in the outer disc, is also contributing to increased metallicity at large galactocentric distances.

In addition to gas flows, stellar migration presents another mechanism that has been shown to flatten stellar metallicity gradients both in N-body simulations and in large cosmological simulations (e.g. Roškar et al. 2008; Di Matteo et al. 2013; Grand, Kawata & Cropper 2015; El-Badry et al. 2016; Buck 2020; Vincenzo & Kobayashi 2020). Nevertheless, the extent to which stellar migration happen in real galaxies has proven to be difficult to assess in observations so far.

Finally, inside-out quenching – either due to morphological quenching (Martig et al. 2009), central compaction (Tacchella et al. 2016), or feedback from active galactic nuclei (Guo et al. 2019) – constitutes another effective way of halting the chemical evolution in the centres of galaxies first while outer parts continue to increase their metal content (e.g. Lacerna et al. 2020).

In summary, there are a variety of processes seen in simulations that affect the chemical evolution of galaxies and are able to substantially alter the radial distribution of metals. Our observations presented in this paper show clear evidence for a primarily mass density-driven metal distribution. This distribution, nevertheless, requires additionally a process (or conglomeration of processes) that promotes enhanced chemical enrichment at larger galactic radii.

In order to allow for a better connection of the physical processes in simulations with observations and to quantitatively compare the predicted and observed distributions of stellar populations in galaxies, it will be interesting to explore whether this observed |$r\Sigma _\star ZR$| can be reproduced in mock observations of cosmological simulations.

6 CONCLUSIONS

This work constitutes a detailed analysis of spatially resolved stellar populations using a very large sample of 7439 galaxies from the MaNGA survey. The sample size allows us to split galaxies by mass and morphology and to push simultaneously the statistical uncertainties to very low values. We focus specifically on the local variations of stellar metallicity with stellar surface mass density and galactocentric distance in the global mass–morphology plane. Our results can be summarized as follows.

• We successfully reproduce the global age–mass and mass–metallicity relations for galaxies of different morphologies.

• We present radial age, metallicity, and surface mass density distributions and find that metallicity gradients are mostly negative but vary strongly with mass, while age profiles are mostly mildly positive in ETGs and negative but shallow in LTGs. Surface mass density profiles are consistently negative with steeper gradients for more massive galaxies.

• We explore in detail the spatially resolved surface mass density–metallicity relation |$r\Sigma _\star ZR$| and find that mass density drives stellar metallicity locally, which connects negative metallicity gradients – as found in most galaxies – to negative surface mass density gradients. In addition, the |$r\Sigma _\star ZR$| is modulated by radial distance.

• In particular, we find that at fixed mass density, metallicity increases with radius. This holds true for the largest part of the sample below log M/M_⊙ ∼ 10.8 except for the lowest mass Es. For log M/M_⊙ > 10.8 the positive trend is only seen in spirals at high surface mass densities, while at low Σ_⋆ the correlation is slightly negative. The positive radial dependence is strongest for spiral galaxies.

• The previous result requires an additional driver of stellar metallicity that promotes chemical enrichment in the outer parts of galaxies more strongly than in the inner parts. We discuss gas accretion, outflows, quenching, and stellar migration as possibilities. A more direct comparison of the effect of these processes as seen in cosmological simulations on the distribution of stellar populations as seen in the observation will help to shed further light on that matter.

ACKNOWLEDGEMENTS

We thank the anonymous referee for comments that helped to improve the paper. The Science and Technology Facilities Council is acknowledged for support through the Consolidated Grant Cosmology and Astrophysics at Portsmouth, ST/S000550/1. JL is supported by the National Science Foundation under Grant No. 2009993. JKB-B acknowledges support from the grant IA-100420 (DGAPA-PAPIIT, UNAM), and funding from the CONACYT grants CF 19-39578, CB-285080, and FC-2016-01-1916.

Numerical computations were done on the Sciama High Performance Compute (HPC) cluster that is supported by the ICG, SEPnet, and the University of Portsmouth.

Funding for the Sloan Digital Sky Survey IV has been provided by the Alfred P. Sloan Foundation, the U.S. DOE Office of Science, and the Participating Institutions. SDSS acknowledges support and resources from the Center for High-Performance Computing at the University of Utah. The SDSS web site is www.sdss.org.

SDSS is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS Collaboration including the Brazilian Participation Group, the Carnegie Institution for Science, Carnegie Mellon University, the Chilean Participation Group, the French Participation Group, Harvard-Smithsonian Centre for Astrophysics, Instituto de Astrofísica de Canarias, The Johns Hopkins University, Kavli Institute for the Physics and Mathematics of the Universe (IPMU)/University of Tokyo, the Korean Participation Group, Lawrence Berkeley National Laboratory, Leibniz Institut für Astrophysik Potsdam (AIP), Max-Planck-Institut für Astronomie (MPIA Heidelberg), Max-Planck-Institut für Astrophysik (MPA Garching), Max-Planck-Institut für Extraterrestrische Physik (MPE), National Astronomical Observatories of China, New Mexico State University, New York University, University of Notre Dame, Observatório Nacional/MCTI, The Ohio State University, Pennsylvania State University, Shanghai Astronomical Observatory, United Kingdom Participation Group, Universidad Nacional Autónoma de México, University of Arizona, University of Colorado Boulder, University of Oxford, University of Portsmouth, University of Utah, University of Virginia, University of Washington, University of Wisconsin, Vanderbilt University, and Yale University.

DATA AVAILABILITY

We shall release the data underlying this paper only through the official releases of the SDSS-IV.

Footnotes

1

There is an overlap between the main sample and the MaNGA ancillary programmes. This first selection only requires the galaxy to be part of the main sample and it does not deselect galaxies that are also part of an ancillary programme.

2

M. Blanton; www.nsatlas.org.

3

https://www.icg.port.ac.uk/firefly/

4

Note that the DR17 version of this catalogue has not yet officially been released. See https://www.sdss.org/dr15/manga/manga-data/manga-firefly-value-added-catalog/ for the latest public release of the MaNGA firefly VAC. A final release is planned for the upcoming SDSS DR17 by the end of 2021.

5

In the literature sometimes also called |$r \rm MZR$|⁠, |$r\mu _\star ZR$|⁠, or Σ_⋆–[Z/H].

APPENDIX: EXAMPLE 2D MAPS FROM THE MaNGA firefly VAC

In Fig. A1, we show example maps of 2D stellar population properties from the MaNGA firefly VAC.

Figure A1.

Example 2D maps from the MaNGA firefly VAC.