Signatures of stellar accretion in MaNGA early-type galaxies

TL;DR

This study uses MaNGA data to analyze stellar metallicity profiles in early-type galaxies, revealing that metallicity flattening in outer regions indicates stellar accretion, with higher ex-situ fractions in more massive galaxies.

Contribution

It demonstrates that radial metallicity profiles are complex and that their flattening in galaxy outskirts serves as evidence for stellar accretion, providing new insights into galaxy assembly.

Findings

Metallicity profiles are non-linear and poorly characterized by a single gradient.

Outer metallicity flattening correlates with stellar accretion signatures.

Ex-situ stellar mass fractions increase with galaxy mass, reaching 80% in the most massive ETGs.

Abstract

The late assembly of massive galaxies is thought to be dominated by stellar accretion in their outskirts (beyond 2 effective radii Re) due to dry, minor galaxy mergers. We use observations of 1010 passive early-type galaxies (ETGs) within z<0.15 from SDSS IV MaNGA to search for evidence of this accretion. The outputs from the stellar population fitting codes FIREFLY, pPXF, and Prospector are compared to control for systematic errors in stellar metallicity (Z) estimation. We find that the average radial logZ/Zsun profiles of ETGs in various stellar mass (M) bins are not linear. As a result, these profiles are poorly characterized by a single gradient value, explaining why weak trends reported in previous work can be difficult to interpret. Instead, we examine the full radial extent of stellar metallicity profiles and find them to flatten in the outskirts of M>10^{11}Msun ETGs. This is a signature of stellar accretion. Based on a toy model for stellar metallicity profiles, we infer the ex-situ stellar mass fraction in ETGs as a function of M and galactocentric radius. We find that ex-situ stars at 2Re make up 20% of the projected stellar mass of M<10^{10.5}Msun ETGs, rising up to 80% for M>10^{11.5}Msun ETGs.

Figures

Tables

Table 1. Table 1: Priors used in our Prospector runs

Parameter	Prior
$\tau$	$\text{LogUniform}({10}^{-2},10)$
dust2	$\text{TopHat}(0,1)$
Stellar age [Gyr]	$\text{TopHat}(5,14)$
Stellar metallicity [logZ/Z⊙]	$\text{TopHat}(-2,0.3)$
Stellar mass [M⊙]	$\text{LogUniform}({10}^5,{10}^{12})$
${\sigma }_{\ast }$ [km/s]	$\text{TopHat}(10,400)$

Table 2. Table 2: Number of spectra used in our analysis a,b a,b footnotemark:

Stellar mass bin [M⊙]	R0.5R${}_{\text{e}}$	0.5R${}_{\text{e}}$R1R${}_{\text{e}}$	1R${}_{\text{e}}$R1.5R${}_{\text{e}}$	1.5R${}_{\text{e}}$R2R${}_{\text{e}}$	2R${}_{\text{e}}$R2.5R${}_{\text{e}}$
${10}^{10}-{10}^{10.5}$	174 (174)	174 (174)	174 (174)	170 (170)	167 (160)
${10}^{10.5}-{10}^{11}$	267 (267)	267 (267)	266 (265)	264 (264)	252 (246)
${10}^{11}-{10}^{11.5}$	420 (420)	420 (420)	417 (417)	393 (392)	347 (319)
${10}^{11.5}-{10}^{12}$	148 (148)	147 (147)	142 (142)	114 (111)	88 (67)

Equations

$$\text{mbox}{logZ}_{\text{mbox}obs}(\text{mbox}{R,lowM}_{\ast })$$

$$\text{mbox}{logZ}_{\text{mbox}in-situ}(\text{mbox}{R,M}_{\ast })=\text{mbox}{logZ}_{\text{mbox}obs}(\text{mbox}{R,lowM}_{\ast })$$

$$\displaystyle+\mbox{logZ}_{\mbox{\scriptsize obs}}(\mbox{0.25R${}_{\mbox{\scriptsize e}}$, M}_{*})-\mbox{logZ}_{\mbox{\scriptsize obs}}(\mbox{0.25R${}_{\mbox{\scriptsize e}}$, low M}_{*})$$

$$\displaystyle\mbox{logZ}_{\mbox{\scriptsize ex-situ}}(\mbox{M}_{*})=\mbox{logZ}_{\mbox{\scriptsize obs}}(\mbox{0.25R${}_{\mbox{\scriptsize e}}$, M}_{*})-\varepsilon$$

$$\text{mbox}{logZ}_{\text{mbox}obs}(\text{mbox}{R,M}_{\ast })=f_{\text{mbox}in-situ}(\text{mbox}{R,M}_{\ast })\text{mbox}{logZ}_{\text{mbox}in-situ}(\text{mbox}{R,M}_{\ast })$$

$$+f_{\text{mbox}ex-situ}(\text{mbox}{R,M}_{\ast })\text{mbox}{logZ}_{\text{mbox}ex-situ}(\text{mbox}{M}_{\ast })$$

$$\text{mbox}[MgFe{]}^{\prime }=\sqrt{\text{mbox}Mgb(0.72\times \text{mbox}Fe5270+0.28\times \text{mbox}Fe5335)}$$

$$\displaystyle\mbox{$<$Fe$>$}=0.5\times\mbox{Fe5270}+0.5\times\mbox{Fe5335}$$

Full text

1 Introduction

The effective radii (R${}_{\mbox{\scriptsize e}}$) of $z\sim 0$ early-type galaxies (ETGs) are observed to be a factor of three to six larger than those of their $z\sim 2$ counterparts (Toft et al. 2007; Cimatti et al. 2008; Buitrago et al. 2008; van Dokkum et al. 2010). On the other hand, the stellar masses (M∗) of local ETGs have only increased by a factor of two since $z\sim 2$ (Daddi et al. 2005; Trujillo et al. 2006a, b, 2007; Zirm et al. 2007; van der Wel et al. 2008; van Dokkum et al. 2008; Damjanov et al. 2009; Cassata et al. 2010, 2011). While galaxies quenched at later times tend to be larger, driving the average R${}_{\mbox{\scriptsize e}}$ upward (progenitor bias; e.g. Valentinuzzi et al. 2010; Carollo et al. 2013), this alone is not sufficient to explain size growth (e.g. Furlong et al. 2017). Late stellar accretion in spheroidal, or even disk configurations (Graham et al. 2015), appears to be required, especially at the high M∗ end (M${}_{*}>$1010.5M*⊙*, e.g. Genel et al. 2018). Minor mergers have been shown to be particularly efficient at increasing the R${}_{\mbox{\scriptsize e}}$ of ETGs while keeping their M∗ roughly constant (e.g. Bezanson et al. 2009; Hopkins et al. 2010; Barro et al. 2013; Cappellari et al. 2013; Wellons et al. 2016).

These ideas are at the basis of the current cosmological picture for structure evolution at $z<2$, in which massive systems accrete stellar envelopes from satellite galaxies (Oser et al. 2010, 2012; Johansson et al. 2012; Moster et al. 2013; Furlong et al. 2017). In this framework, stars that formed within their host galaxies tend to dominate at the center, whereas accreted stars begin to do so in the outskirts (R$\sim$2R${}_{\mbox{\scriptsize e}}$; Rodriguez-Gomez et al. 2016) and in the lower surface brightness regions beyond 2R${}_{\mbox{\scriptsize e}}$ known as stellar halos (Zolotov et al. 2009; Tissera et al. 2013, 2014; Cooper et al. 2015). These stellar populations of different origin are usually referred to as in-situ and ex-situ, respectively. Several simulations have made predictions about observational signatures of the predicted radial transition from in-situ to ex-situ (e.g. Pillepich et al. 2014). Among stellar population tracers, stellar metallicity is expected to be one of the most sensitive to this transition (e.g. Cook et al. 2016).

In the absence of late-time minor mergers, the radial stellar metallicity profiles are predicted to be negative (Kobayashi 2004; Pipino et al. 2010; Taylor & Kobayashi 2017). This implies that the outer parts of ETGs tend to be more metal-poor than the inner parts. Albeit with significant variance (Lackner et al. 2012; Hirschmann et al. 2015), the deposition of accreted stars in the outskirts of galaxies induces flattening of the in-situ profile (Cook et al. 2016; Taylor & Kobayashi 2017). Since mergers are expected to have a larger effect on more massive systems, the resulting prediction is that the stellar metallicity profiles of ETGs are flatter toward higher M∗, especially in the stellar halos (Cook et al. 2016).

These theoretical predictions have motivated the search for observational signatures of stellar accretion. Using long-slit spectroscopy, Carollo et al. 1993 estimated the strength of metallicity-sensitive stellar absorption features as a function of galactocentric radius in 42 nearby galaxies. Though larger samples can be studied using photometric surveys (e.g. La Barbera et al. 2005, 2011; Tortora et al. 2010; Tortora & Napolitano 2012), spectroscopy is critical for breaking the age-metallicity degeneracy. More recently, studies of stellar populations in nearby galaxies have benefited from integral field unit (IFU) surveys like MASSIVE (Greene et al. 2013, 2015), CALIFA (Sánchez et al. 2012), SAMI (Allen et al. 2015), and MaNGA (Mapping Nearby Galaxies at Apache Point Observatory; Bundy et al. 2015). In particular, MaNGA observations extend to the outskirts of galaxies (beyond 2 R${}_{\mbox{\scriptsize e}}$), starting to probe the radii at which the signatures of minor mergers are predicted to appear (e.g. Cook et al. 2016).

Stellar metallicity profiles are typically characterized by radial gradients, estimated by fitting a linear form to the profile between the center and 1-2R${}_{\mbox{\scriptsize e}}$ (e.g. Zheng et al. 2017; Goddard et al. 2017b; Li et al. 2018). In agreement with simulations, the metallicity gradients of ETGs tend to be negative (e.g. Rawle et al. 2010; González Delgado et al. 2015; Roig et al. 2015). However, the dependence of the gradient slope on M∗ remains unclear. Based on a sample of $\sim 10^{3}$ galaxies from the MaNGA survey, Zheng et al. (2017) find weak or no correlation between the gradients and M∗. Using data from the same survey, Goddard et al. (2017a) find that gradients are steeper with increasing M∗, although with low significance. Though also based on MaNGA, Li et al. (2018) find shallower gradients at higher central velocity dispersions ($\sigma_{*}>100$km/s). There are several possible sources for these discrepancies, from stellar population synthesis approach (see Conroy 2013) to fitting method. Another important factor, as we show in this paper, is that the stellar metallicity profiles of ETGs are not well described by a linear fit.

In this work, we examine the full radial extent of metallicity profiles from spatially resolved spectroscopy of 1010 ETGs from MaNGA. We inform our interpretation of the stellar metallicity profiles by using results from hydrodynamical simulations (e.g. Cook et al. 2016; Rodriguez-Gomez et al. 2016; D’Souza & Bell 2018). This paper is structured as follows. In Section 2 we define our sample. In Section 3, we describe the stellar population fitting process with the codes FIREFLY (Wilkinson et al. 2017), pPXF (Cappellari & Emsellem 2004; Cappellari 2017), and Prospector (Leja et al. 2017). We show our results in Section 4 and discuss the implications in Section 5. We summarize in Section 6. This work adopts H${}_{0}=70$ km s*-1Mpc-1* and all magnitudes are reported in the AB system (Oke & Gunn 1983).

2 Dataset

The MaNGA survey (Bundy et al. 2015; Yan et al. 2016a) is part of the fourth generation of SDSS (York et al. 2000; Gunn et al. 2006; Blanton et al. 2017), and is on track to provide spatially resolved spectra for ten thousand nearby galaxies ($z<0.15$) by the end of 2020. By means of integral field unit spectroscopy (IFS; Smee et al. 2013; Drory et al. 2015; Law et al. 2015), every galaxy is observed with 19-to-127 fiber bundles with diameters varying between 12$\farcs$5 and 32$\farcs$5. The resulting radial coverage reaches between 1.5R${}_{\mbox{\scriptsize e}}$ and 2.5R${}_{\mbox{\scriptsize e}}$ for most targets (Wake et al. 2017; see Figure 1). The spectra cover the wavelength range 3600-10300 Å at a resolution of R$\sim$2000.

All MaNGA data used in this work were reduced by the Data Reduction Pipeline (DRP; Law et al. 2016; Yan et al. 2016b). The reduced spectra have a median spectral resolution of $\sigma$=72 km s*-1*. The data cubes typically reach a 10$\sigma$ continuum surface brightness of $\mu$=23.5 mag arcsec*-2*, and their astrometry is measured to be accurate to 0$\farcs$1 (Law et al. 2016). De-projected distances and stellar kinematic maps have been calculated by the MaNGA Data Analysis Pipeline (DAP; Westfall et al. 2019). This work also makes use of Marvin (Cherinka et al. 2017), the specially designed tool for access and handling of MaNGA data111https://api.sdss.org/doc/manga/marvin.

This paper is based on the SDSS Data Release 15 (DR15), which consists of the observations of the first 4675 MaNGA targets. We extracted the stellar masses (M∗), Sersic indices ($n_{\mbox{\scriptsize Sersic}}$), and effective radii (R${}_{\mbox{\scriptsize e}}$) of these galaxies from the publicly available NASA-Sloan Atlas222http://nsatlas.org(NSA). In particular, the M∗ estimates were derived using a k-correction fit to the Sersic fluxes (Blanton & Roweis 2007), adopting the Bruzual & Charlot (2003) stellar population models and a Chabrier (2003) initial mass function (IMF). They also assumed H${}_{0}=100$ km s*-1Mpc-1*, but we scaled them for an H${}_{0}=70$ km s*-1Mpc-1* cosmology. The $n_{\mbox{\scriptsize Sersic}}$ estimates were obtained from one-component, two-dimensional fits to r-band images. The R${}_{\mbox{\scriptsize e}}$ are determined using an elliptical Petrosian analysis of the $r$-band image from the NSA. All NSA measurements use the detection and deblending technique described in Blanton et al. (2011).

To select ETGs, we first applied the morphological cut $n_{\mbox{\scriptsize Sersic}}>2.5$ (e.g. Blanton et al. 2003, 2005; Peng et al. 2010). In addition, we selected passive ETGs by using the average H$\alpha$ equivalent width across the galaxy -EW(H$\alpha)$- as proxy for specific star-formation rate (sSFR). The cut was EW(H$\alpha)<3$ Å, which is commonly used to distinguish between ionization due to smooth background of hot evolved stars and due to star formation and AGN (Cid Fernandes et al. 2011; see also Belfiore et al. 2016). This yielded a sample with 1101 galaxies. We also limited the central velocity dispersions and stellar masses of our sample to the ranges $\sigma_{*}<400$ km s*-1* and $10<$logM/M${}_{*}<12$, respectively. We performed these cuts to provide a relatively uniform distribution of ETGs over M∗. The final outcome was a sample of 1010 ETGs. We did not remove quiescent galaxies with significant stellar disks from the sample. From visual inspection, we estimate the fraction of lenticulars (S0s) to be $\lesssim 20$ %. However, we acknowledge the challenge of achieving precise S0 classification of SDSS galaxies (see Nair & Abraham 2010). Our selection may also miss blue ellipticals, but their number fraction is $\lesssim 5\%$ for our M∗ range (Kannappan et al. 2009). Our goal here is to study a generally passive sample of spheroidal galaxies. We delay to future work a characterization of stellar populations in more finely discriminated morphological types.

3 Methodology

3.1 Radial binning

Using the R${}_{\mbox{\scriptsize e}}$ value of every galaxy, we associated elliptical polar radii to all spaxels in units of R${}_{\mbox{\scriptsize e}}$. These account for the axis ratio of every object, which were measured on the r-band photometry. We then binned them into the five annuli R/R${}_{\mbox{\scriptsize e}}$= $[0,0.5]$, $[0.5,1]$, $[1,1.5]$, $[1.5,2]$, and $[2,2.5]$. This is shown for a sample galaxy on the left panel of Figure 1.

After binning, we shifted every spectrum back to the rest-frame using the stellar systemic velocity ($v_{*}$) maps calculated by the DAP. We used the maps computed with a Voronoi binning scheme that aims for a minimum signal-to-noise ratio of 10 per bin. For each galaxy, we co-added the spectra in every annular bin. We did not convolve the spectra to a common $\sigma_{*}$ prior to stacking. After co-addition, we ran pPXF (Cappellari & Emsellem 2004; Cappellari 2017) with the MILES Single Stellar Population (SSP) library (Vazdekis et al. 2010) on the stacked spectra to measure the co-added $v_{*}$ and $\sigma_{*}$. The right panel of Figure 1 shows the five co-added spectra for a sample galaxy.

3.2 Stellar population fitting

Estimates of stellar population parameters like stellar metallicity can be obtained by full spectral fitting, but depend sensitively on the adopted priors, assumptions used to generate template spectra (Conroy 2013), and fitting method. To mitigate the effect of systematic biases from any one approach, we applied three independent codes to the same data and examine the differences that arise.

The first code we ran was the public version of FIREFLY333[FIREFLY - A full spectral fitting code

http://www.icg.port.ac.uk/FIREFLY/](http://www.icg.port.ac.uk/FIREFLY/)444[https://github.com/FireflySpectra/firefly_release](https://github.com/FireflySpectra/firefly_release)(Comparat et al. 2017; Wilkinson et al. 2017; Goddard et al. 2017a). This $\chi^{2}$ minimization code decouples stellar populations from dust by removing the low-order continuum shape before performing the model fitting. Hence, it focuses on high frequency modes in the spectra to infer stellar ages and metallicities. SSPs of different ages and metallicities are added iteratively until the improvement in $\chi^{2}$ is negligible.

We ran the code with the stellar population models of Maraston & Strömbäck (2011), MILES stellar library (Sánchez-Blázquez et al. 2006), and Chabrier IMF (Chabrier 2003). We used a set of SSPs covering an age grid between 6.5 Myr and 15 Gyr, while the sampled stellar metallicities were logZ/Z⊙=-2.3, -1.3, -0.3, 0.0, and 0.3. The library spans the wavelength range 4000Å to 7400Å. As shown in Wilkinson et al. 2017, FIREFLY effectively recovers stellar population parameters for spectra with S/N$>10$ (see also Goddard et al. 2017b). To limit the systematics in the measurements from Firefly, we excluded any co-added spectra with S/N$<10$. We also masked emission lines. Fitting with FIREFLY took, on average, a minute per spectrum on a single core. Throughout this paper, we show light-weighted measurements, although we find similar results when using the mass-weighted counterparts.

We also ran pPXF555[pPXF

https://www-astro.physics.ox.ac.uk/ mxc/software/#ppxf](https://www-astro.physics.ox.ac.uk/%C2%A0mxc/software/#ppxf)(Cappellari & Emsellem 2004; Cappellari 2017) on our spectra. This code applies a penalized maximum likelihood approach to fit libraries of stellar population templates to observed data. Since this code penalizes pixels that are not well characterized by the templates, it minimizes template mismatch. We ran it with the included library of SSPs based on the MILES stellar library (Sánchez-Blázquez et al. 2006; Vazdekis et al. 2010).

We simultaneously fitted for the gas and the stars, allowing for two moments in gas kinematics and four in stellar kinematics. We chose not to smooth the distribution of template weights (i.e., no regularization). After the best linear combination of templates was found, we added several realizations of the noise in the spectra to the best fit. This allowed us to characterize the uncertainties in the reported stellar population parameters. On average, our runs of pPXF took about a minute per spectrum on a single core.

The third stellar population fitting code we ran was Prospector666[Prospector

https://github.com/bd-j/prospector/blob/master/doc/index.rst](https://github.com/bd-j/prospector/blob/master/doc/index.rst)(Leja et al. 2017). This code is based on the stellar population synthesis code FSPS777[FSPS: Flexible Stellar Population Synthesis

https://github.com/cconroy20/fsps](https://github.com/cconroy20/fsps)(Conroy et al. 2009; Conroy & Gunn 2010), which generates composite stellar spectra for a variety of prescriptions for stellar population synthesis and evolution. This allows Prospector to sample the posterior distribution of a user-defined parameter space, while formally characterizing uncertainties and degeneracies. We chose the MILES stellar library (Sánchez-Blázquez et al. 2006), MIST isochrones (Dotter 2016; Choi et al. 2016), and Kroupa IMF (Kroupa 2001) as inputs. We also masked emission lines prior to fitting.

Since we fitted old stellar populations, we modeled the spectra with exponentially decaying ($\tau$) star-formation histories to speed up the fitting process. In addition to $\tau$, our parameter space included the optical depth of dust in the V-band, and stellar ages, metallicities, masses, and velocity dispersions. Our priors are shown in Table 1. To derive the posterior distributions, we used the Dynamic Nested Sampling package dynesty888[dynesty

https://github.com/joshspeagle/dynesty/blob/master/docs

/source/index.rst](https://github.com/joshspeagle/dynesty/blob/master/docs/source/index.rst)(Speagle 2019). On average, convergence of Prospector with dynesty was achieved after an hour per spectrum on a single core.

4 Results

Using the three codes described above, we derived stellar population parameters in each radial bin for all galaxies in the sample. After binning in M∗ (with numbers in Table 2), we computed the average stellar metallicity profiles as a function of M∗ and show them in Figure 2. The three panels show the results from the three fitting codes. While the metallicity profiles differ in normalization and in their detailed shapes, qualitative trends are similar across the codes.

We start by discussing the two notable discrepancies among the outputs. First, pPXF systematically measures metallicities $\sim 0.1$ dex lower than Firefly and Prospector. This overall offset does not correlate with S/N or M∗ and will not affect our primary conclusions, which are based on the shape of derived metallicity profiles. Second, Firefly outputs tend to avoid metallicities in the range logZ/Z${}_{\odot}=[-1.3,-0.3]$, preferring higher values. This is presumably due to sampling in the stellar metallicity grid (see Wilkinson et al. 2017). As we show in Figure 2, flattening of the Firefly metallicity profiles occurs at higher metallicities as a result.

In nearly all radial bins, more massive galaxies exhibit more metal-rich stars. The logZ/Z*⊙* profiles of ETGs fall linearly with galactocentric radius out to 1.5R${}_{\mbox{\scriptsize e}}$. Remarkably, the profiles flatten at the largest radii for M${}_{*}>10^{11}$M*⊙*. The flattening is present in the output of all three codes. Comparing a given set of profiles as a function of M∗, we see that the radius at which this flattening occurs moves inward as M∗ increases. These results are also apparent in the behavior of Lick indices Fe4531, Mgb(5178), and Fe5270 (Appendix A). The observed flattening is consistent with the signatures of stellar accretion predicted by hydrodynamical simulations (e.g. Cook et al. 2016) and motivates the interpretative framework we discuss in Section 5.

We note that even though the M∗ dependence of the stellar metallicity profiles is consistent across codes, the same cannot be said about the stellar age profiles (not shown). This is not surprising, since it is extremely difficult to determine the ages of stellar populations older than 9 Gyr because of the slow isochrone evolution at late times (Conroy 2013). Since radial gradients in stellar age are not predicted to capture much information about the accretion history of ETGs (Cook et al. 2016), we leave a more detailed analysis of stellar ages for future work.

Some galaxies only satisfied our quality criteria (see Section 3) at some annuli. Hence, some galaxies contributed only to some regions in the profiles of Figure 2. To ensure our results are not biased, as a result we constructed a subset of 822 ETGs composed only of high quality spectra (S/N$>10$ for all radii). Our results were also recovered with this subset.

We have also attempted to reproduce our results using the publicly available Firefly999MaNGA FIREFLY Value Added Catalog http://www.sdss.org/dr14/manga/manga-data/manga- FIREFLY-value-added-catalog/ (Goddard et al. 2017a, b) and Pipe3D101010Pipe3D Value Added Catalog https://www.sdss.org/dr14/manga/manga-data/manga-pipe3d-value-added-catalog/ (Sánchez et al. 2016, 2018) Value Added Catalogs, which provide spatially-resolved maps of stellar population properties for MaNGA galaxies. Unfortunately, Voronoi bins with S/N$<$10 dominate in the outermost low-surface brightness regions. Various tests have shown that stellar population codes are biased at S/N$<$10 (e.g. Wilkinson et al. 2017). As a result of these complications, we refrained from incorporating these catalogs in our analysis.

5 Discussion

5.1 On the radial metallicity profiles of ETGs

In the R$<$R${}_{\mbox{\scriptsize e}}$ region, Martín-Navarro et al. (2018) found that the stellar metallicity profiles of ETGs fall more steeply at higher $\sigma_{*}$ and M∗. Similarly, Goddard et al. (2017a) reported weak evidence for a steepening of their radial gradients with M∗. On the other hand, Kuntschner et al. (2010); Tortora et al. (2010); Kuntschner (2015); Li et al. (2018) found gradients to flatten at higher $\sigma_{*}$. González Delgado et al. (2015); Zheng et al. (2017) claimed no clear correlation between their stellar metallicity gradients and M∗. Similarly, Greene et al. (2013, 2015) found no strong correlations between the shape of element abundance profiles and $\sigma_{*}$. In this work, we found the profiles to flatten in the outskirts for logM∗/M${}_{\odot}\gtrsim 11$. Here, we demonstrate how some of the apparent disagreement among observations may owe to the definition of metallicity gradients.

A quick look at our Figure 2 reveals that the average metallicity profiles of ETGs are not straight lines. It stands to reason that fitting lines to these radial profiles could “wash-out” the flattening in the outskirts of high M∗ ETGs. Figure 3 shows the outcome of fitting lines to our metallicity profiles over different radial ranges motivated by the literature. Some ranges trace the inner regions (R$<$Re; Li et al. 2018), while others have more extended coverage (R$<$2Re; Goddard et al. 2017a). The scatter is considerable in all cases, and recovering any correlations with M∗ is difficult. We conclude that gradients are sensitive to radial coverage (see also Greene et al. 2019) and can also miss important behavior in the stellar metallicity profiles. Gradients should be avoided when possible.

5.2 Comparison with hydrodynamical simulations

Hydrodynamical simulations predict stellar accretion to induce gradient flattening (e.g. Cook et al. 2016). In general, stars accreted via dry, minor mergers tend to settle around and beyond the outskirts of ETGs (R=2-4R${}_{\mbox{\scriptsize e}}$), which results in a flatter stellar metallicity profile than the inherently steeper form it originally had. Since mergers are expected to have a larger effect on more massive systems, this prediction is in broad agreement with our results from Figure 2.

A relevant point involves the radii at which accretion signatures are expected to appear. Rodriguez-Gomez et al. (2016) derived the accreted mass fraction of galaxies as a function of galactocentric radius in the Illustris simulation (Vogelsberger et al. 2014a, b). On average, this fraction increases with radius. It goes from zero at the center to unity at radii R$\gtrsim 5$R${}_{\mbox{\scriptsize e}}$. This motivates the definition of the transition radius (R${}_{\mbox{\scriptsize T}}$). It is defined as the galactocentric radius at which the M∗ fraction of the ex-situ stellar component overtakes its in-situ counterpart (D’Souza et al. 2014). Rodriguez-Gomez et al. (2016) found R${}_{\mbox{\scriptsize T}}$ to decrease with M∗, going from R${}_{\mbox{\scriptsize T}}\sim 5$R${}_{\mbox{\scriptsize e}}$ at M${}_{*}\sim 10^{10}$M*⊙* to R${}_{\mbox{\scriptsize T}}<$R${}_{\mbox{\scriptsize e}}$ at M${}_{*}\sim 10^{12}$M*⊙*. Our results are qualitatively consistent with this prediction.

However, there are some quantitative tensions. For logM∗/M${}_{\odot}\sim 11$ galaxies, Rodriguez-Gomez et al. (2016) reported R${}_{\mbox{\scriptsize T}}\sim 4$ R${}_{\mbox{\scriptsize e}}$. Within 2.5R${}_{\mbox{\scriptsize e}}$, we should only be probing accreted stellar mass fractions of $\lesssim$0.3 at this mass range. Cook et al. (2016), also based on the Illustris simulation, reported that the flattening of metallicity gradients with M∗ only becomes noticeable in the stellar halo (R$=2-4$R${}_{\mbox{\scriptsize e}}$). Therefore, the signatures we see in Figure 2 are apparent at smaller radii than some simulations have predicted. There are a few possible explanations for this tension. The works of Rodriguez-Gomez et al. (2016); Cook et al. (2016) were based on Illustris. Galaxies at $z\sim 0$ from the first generation of this simulation were found to be larger by a factor of $\sim 2$ than observed galaxies. IllustrisTNG solved this problem, among others, by improving the treatment of galactic winds, magnetic fields, and black hole feedback (Pillepich et al. 2018; Weinberger et al. 2017, 2018). The treatment of these, among with other secular processes, can strongly impact the stellar population gradients measured in simulations (e.g. Taylor & Kobayashi 2017). On the observational side, estimates of the ages and metallicities of stellar populations can strongly depend on the choice of stellar library, isochrones, and approach to fitting. These systematic uncertainties also affect the conversion between stellar mass and stellar light, impacting the comparison between simulations and observations.

5.3 Estimating the ex-situ stellar mass fraction

Observationally, global stellar metallicity correlates with M∗ or the central velocity dispersion $\sigma_{*}$ of galaxies (Faber & Jackson 1976; Cid Fernandes et al. 2005; Gallazzi et al. 2005; Thomas et al. 2005, 2010; González Delgado et al. 2014), as would be expected if the deeper potential wells of more massive systems limit the impact of galactic winds (Matteucci 1994). In the Illustris simulation, D’Souza & Bell (2018) found an accreted M${}_{\mbox{\scriptsize acc}}$-Z${}_{\mbox{\scriptsize acc}}$ relation, where M${}_{\mbox{\scriptsize acc}}$ and Z${}_{\mbox{\scriptsize acc}}$ refer to the stellar mass and stellar metallicity of the accreted components, respectively. This relationship lies $\sim 0.3$ dex below the global counterpart. We can make informed assumptions for the in-situ stellar metallicity profile and the M${}_{\mbox{\scriptsize acc}}$-Z${}_{\mbox{\scriptsize acc}}$ relation to build a toy model capable of inferring the ex-situ M∗ fraction as a function of mass and galactocentric radius from our observations.

We assume the intrinsic in-situ metallicity profiles of ETGs to be well described by the profiles observed in the low mass end of our sample. This is supported by hydrodynamical simulations that find M${}_{*}\sim 10^{10}$M*⊙* galaxies to be dominated by in-situ stars within the radial coverage of our data (Rodriguez-Gomez et al. 2016). We take the M${}_{*}=10^{10}-10^{10.5}$M*⊙* profiles from Figure 2 for each code and refer to them as

[TABLE]

In our model, the in-situ profiles of all galaxies follow the shape of $\mbox{logZ}_{\mbox{\scriptsize obs}}(\mbox{R, low M}_{*})$ with a normalization applied to match the metallicity at the center (i.e. within 0.5 R${}_{\mbox{\scriptsize e}}$). This can be written as

[TABLE]

with a schematic representation in Figure 4.

The M${}_{\mbox{\scriptsize acc}}$-Z${}_{\mbox{\scriptsize acc}}$ relation is offset 0.3 dex from the global counterpart in the Illustris simulation. The existence of this relation originates from single massive progenitors contributing to the bulk of the mass to the accreted stellar component (D’Souza & Bell 2018). If we assume the accreted envelopes of ETGs to be comparable in stellar mass to their host ETG (e.g. Rodriguez-Gomez et al. 2016), ex-situ metallicities can be approximated by

[TABLE]

i.e., stellar metallicity of ex-situ stars will be $\varepsilon=0.3$ dex lower than the metallicity at the center of the galaxy. Note that there is no dependence on galactocentric radius in the definition of $\mbox{logZ}_{\mbox{\scriptsize ex}}$.

For measured metallicities, the offset will be dependent on the stellar population synthesis approach. To account for differences between codes, we set $\varepsilon$ equal to the difference in metallicity between the centers and outskirts of M${}_{*}=10^{11.5}-10^{12}$M*⊙* ETGs (see Figure 4). The corresponding values are $\varepsilon\sim-0.14$ (Firefly), $-0.29$ (pPXF), and $-0.24$ (Prospector).

We can now write observed metallicities as a linear combination between in-situ and ex-situ metallicities:

[TABLE]

where $f_{\mbox{\scriptsize in-situ}}$ and $f_{\mbox{\scriptsize ex-situ}}=1-f_{\mbox{\scriptsize in-situ}}$ are the in-situ and ex-situ fractions. Figure 4 describes our toy model and how we derive ex-situ fractions from it.

The results as a function of M∗, galactocentric radius, and code are shown in Figure 5. Ex-situ fractions increasingly dominate at larger radii and higher M∗. Ex-situ stars at R$\sim$2R${}_{\mbox{\scriptsize e}}$ make up $\lesssim$20% of the projected stellar mass of M${}_{*}\lesssim 10^{10.5}$M*⊙* ETGs, rising up to $\gtrsim$80% for M${}_{*}\gtrsim 10^{11.5}$M*⊙* ETGs.

Stellar accretion and minor mergers provide an explanation for the size growth of spheroids from $z\sim 2$ to the present. Keeping in mind the simple nature of our comparison, we showed that the logZ profiles of nearby ETGs are consistent with this framework. However, this picture might not apply to S0s, which we visually estimate to compose $\lesssim 20$% of our sample. The growth and accretion histories of S0s can differ from those of elliptical galaxies ($n_{\mbox{\scriptsize Sersic}}>2.5$; Blanton et al. 2003, 2005; Peng et al. 2010), as suggested by Johnston et al. (2012, 2014); Fraser-McKelvie et al. (2018); Saha & Cortesi (2018) (see also Diaz et al. 2018). Moreover, galaxy assembly history is not only expected to depend on the total M∗ or morphology of galaxies, but also on their environment (e.g. Greene et al. 2015, 2019). We will study second order trends in the metallicity profiles of ETGs in follow-up work.

6 Summary

We characterized the radial stellar metallicity profiles of MaNGA ETGs and compared them with predictions from hierarchical formation. Through stellar population fitting with Firefly, pPXF, and Prospector, we found the following:

1. The three codes are built around different stellar population synthesis codes and are unique in their approach to fitting. Nonetheless, we found the main conclusions from this paper not to be dependent on the fitting code.

2. The profiles of logM∗/M${}_{\odot}\gtrsim 11$ ETGs fall with galactocentric radius and flatten beyond R$\sim$1.5R${}_{\mbox{\scriptsize e}}$. Based on hydrodynamical simulations, a possible explanation for this flattening is stellar accretion through minor mergers.

3. The average radial metallicity profiles of ETGs are not linear. Therefore, linear fits can miss important behavior in the stellar metallicity profiles. When possible, fitting stellar population gradients should be avoided.

4. Using informed assumptions for the in-situ metallicity profile and the metallicity of accreted stars, we built a toy model to infer the ex-situ stellar mass fraction of ETGs. We found ex-situ signatures to grow in significance toward large galactocentric radii and higher M∗.

We thank the referee for their constructive comments and suggestions. We also thank Ben Johnson for helping us run Prospector on MaNGA data. We are grateful to Platon Karpov and Enrico Ramírez-Ruiz for allowing us to use their cluster, Comrade, to run stellar the stellar population fitting codes. This work made use of GNU Parallel (Tange 2018) for executing jobs on Comrade. We thank Brian Cherinka and José Sánchez-Gallego for helping the authors to familiarize themselves with Marvin, the core Python package and web framework for MaNGA data. GO acknowledges support from the Regents’ Fellowship from the University of California, Santa Cruz. KB is supported by the UC-MEXUS-CONACYT Grant. ZZ is supported by the National Natural Science Foundation of China No. 11703036. This research made use of Marvin, a core Python package and web framework for MaNGA data, developed by Brian Cherinka, José Sánchez-Gallego, and Brett Andrews (MaNGA Collaboration, 2018). Funding for the Sloan Digital Sky Survey IV has been provided by the Alfred P. Sloan Foundation, the U.S. Department of Energy 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 Center 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, 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.

Appendix A Lick index profiles

Lick indices (Worthey et al. 1994; Thomas et al. 2003; Parikh et al. 2019) are a useful method to empirically estimate the chemical abundance patterns of galaxies. Here, we compute the radial profiles of Mgb(5178), Fe5270, and Fe5335 to test the high M∗ flattening we find through stellar population fitting. We retrieved the indices measured by the MaNGA Data Analysis Pipeline and used them to compute:

[TABLE]

Here [MgFe]’ and $<$Fe$>$ are tracers of the global and Iron abundances (Johnston et al. 2018). We binned the measurements into the five annuli R/R${}_{\mbox{\scriptsize e}}$= $[0,0.5]$, $[0.5,1]$, $[1,1.5]$, $[1.5,2]$, and $[2,2.5]$ to derive the median profiles shown in Figure 5. Note how the profiles flatten for the highest M∗ bin.