Unveiling the time evolution of chemical abundances across the Milky Way disc with APOGEE

ABSTRACT

Chemical abundances are an essential tool in untangling the Milky Way’s enrichment history. However, the evolution of the interstellar medium abundance gradient with cosmic time is lost as a result of radial mixing processes. For the first time, we quantify the evolution of many observational abundances across the Galactic disc as a function of lookback time and birth radius, |$\rm \text{R}_\text{birth}$|⁠. Using an empirical approach, we derive |$\rm \text{R}_\text{birth}$| estimates for 145 447 APOGEE DR17 red giant disc stars, based solely on their ages and |$\rm [Fe/H]$|⁠. We explore the detailed evolution of six abundances [Mg, Ca (α), Mn (iron-peak), Al, C (light), Ce (s-process)] across the Milky Way disc using 87 426 APOGEE DR17 red giant stars. We discover that the interstellar medium had three fluctuations in the metallicity gradient ∼9, ∼6, and ∼4 Gyr ago. The first coincides with the end of high-α sequence formation around the time of the Gaia–Sausage–Enceladus disruption, while the others are likely related to passages of the Sagittarius dwarf galaxy. A clear distinction is found between present-day observed radial gradients with age and the evolution with lookback time for both [X/Fe] and [X/H], resulting from the significant flattening and inversion in old populations due to radial migration. We find the |$\rm [Fe/H]$|–|$\rm [\alpha /Fe]$| bimodality is also seen as a separation in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| plane for the light and α-elements. Our results recover the chemical enrichment of the Galactic disc over the past 12 Gyr, providing tight constraints on Galactic disc chemical evolution models.

Galaxy: abundances, Galaxy: disc, Galaxy: evolution

1 INTRODUCTION

The precise individual element abundance measurements for large numbers of stars measured by large spectroscopic surveys have advanced the field of Galactic archaeology over the past decade. In particular, the 17th data release of the Apache Point Observatory Galactic Evolution Experiment (APOGEE; Majewski et al. 2017; Abdurro’uf et al. 2022) provides 20+ elemental abundance species for over 650 000 stars. Due to its deep coverage in the near infrared, APOGEE observes stars throughout the Galactic disc (Zasowski et al. 2017), and thus is able to provide a detailed chemical map.

The relatively unchanging¹ chemical abundances of stellar atmospheres give insight into the star’s birth properties (Tinsley 1979; Freeman & Bland-Hawthorn 2002; Minchev et al. 2018; Ratcliffe et al. 2022) and can reveal important evolutionary events during the Milky Way’s lifetime (e.g. Pagel 2009; Matteucci 2012; Frankel et al. 2018; Helmi et al. 2018; Belokurov & Kravtsov 2022; Conroy et al. 2022; Horta et al. 2023; Khoperskov et al. 2022b; Lu et al. 2022a). However, despite the increased number of chemical abundances spanning multiple nucleosynthetic families measured for many stars across the Galactic disc, there are still many unresolved questions. For instance, the origin of the well-known bimodality in the |$\rm [Fe/H]$|–|$\rm [\alpha /Fe]$| plane is still under debate (Grisoni et al. 2017; Mackereth et al. 2018; Clarke et al. 2019; Buck 2020; Lian et al. 2020; Agertz et al. 2021; Khoperskov et al. 2021; Sharma, Hayden & Bland-Hawthorn 2021).

It has also been found that some elements capture similar information, meaning the chemical abundance space collapses on to a lower dimensional manifold which roughly corresponds to the different nucleosynthetic families (Ting et al. 2012; Price-Jones & Bovy 2018; Ness et al. 2019; Weinberg et al. 2019; Ratcliffe et al. 2020; Griffith et al. 2021; Ness et al. 2022). While a lower dimensional abundance space may capture the broad strokes of chemical enrichment, each element uniquely contributes to the chemical evolution of the Galaxy (e.g. Kobayashi, Karakas & Lugaro 2020; Matteucci 2021). In order to best understand the Milky Way’s intricate evolution, it is crucial to describe the way abundances vary across the disc with time. Classical chemical evolution modelling and more modern hydrodynamical simulations have had great success characterizing the Galactic abundance evolution (e.g. Steinmetz & Mueller 1994; Chiappini, Matteucci & Gratton 1997; Minchev, Chiappini & Martig 2013; Kubryk, Prantzos & Athanassoula 2015; Sanders & Binney 2015; Buck 2020; Sanderson et al. 2020; Johnson et al. 2021; Bellardini et al. 2022); however, there is still tension between observed and expected data modelled with current nucleosynthetic yield tables (Rybizki, Just & Rix 2017; Blancato et al. 2019; Buck et al. 2021; Palla et al. 2022). Previous work has utilized open clusters to trace chemical evolution by looking at radial gradients conditioned on age (e.g. Janes 1979; Friel 1995; Chen, Hou & Wang 2003; Magrini et al. 2009, 2023; Netopil et al. 2016; Zhang, Chen & Zhao 2021; Myers et al. 2022; Spina, Magrini & Cunha 2022). However, differences in cluster distance estimates lead to variations in the radial gradient measurements of up to 0.03 dex kpc⁻¹, even when keeping the set of abundance measurements constant (Donor et al. 2018). Individual element abundance gradients have also been measured for field stars using present-day radii (e.g. Mayor 1976; Grenon 1987; Daflon & Cunha 2004; Genovali et al. 2014; Bragança et al. 2019; Eilers et al. 2022) or guiding-centre radii (e.g. Casagrande et al. 2011; Anders et al. 2014).

It is now well established that stars migrate away from their birth sites (e.g. Sellwood & Binney 2002; Roškar et al. 2008; Quillen et al. 2009; Schönrich & Binney 2009a; Minchev & Famaey 2010; Grand, Kawata & Cropper 2012; Frankel et al. 2018; Carr et al. 2022), causing the measured present-day abundance gradients to be flatter than the environment they were born in (Minchev et al. 2012, 2013; Pilkington et al. 2012; Kubryk, Prantzos & Athanassoula 2013; see also fig. 11 in Vincenzo & Kobayashi 2020). Thus, radial migration has begun to wash away chemical signatures in the disc. Knowledge of birth radius (⁠|$\rm \text{R}_\text{birth}$|⁠) can help answer fundamental questions still not yet understood, such as the now-lost evolution of the interstellar medium abundance gradient with cosmic time. With knowledge of birth radii, we would be able to go beyond estimating radial abundance gradients from conditioning on age using present Galactic radii – and therefore getting an incomplete story of the Milky Way’s evolution – and have a clearer understanding of the assembly history of the Galaxy. The origin of the high- and low-α sequences could be understood, in addition to the relative strength of different nucleosynthetic events at specific moments and places in the Galactic disc.

One way to recover birth radii is to leverage the stars’ chemical abundances to recover their birth groups, which are believed to be chemically homogeneous (Bland-Hawthorn, Krumholz & Freeman 2010). However, it has become clear in recent years that the hopes for the feasibility of strong chemical tagging (the idea of directly relating the individual members of dispersed star clusters to their birth sites based only on spectroscopy, in the terminology of Freeman & Bland-Hawthorn 2002) are slim, at least with regards to the Galactic disc. While some studies used the strong chemical tagging technique to find possibly dispersed star clusters (e.g. Price-Jones et al. 2020), it has been shown that even in the most favourable conditions (homogeneous high-resolution spectral analysis of only red clump stars with no contamination from field stars), more than 70  per cent of the recovered statistical groups (using the most advanced clustering methods) contain stars belonging to different real clusters (Casamiquela et al. 2021). While strong chemical tagging is thus discarded for the time being, the prospects for weak chemical tagging in the disc are extremely favourable due to the ever-growing amount of data from high-resolution spectroscopic surveys.

Minchev et al. (2018) devised a method for obtaining stellar birth radii in any observational sample with good age and metallicity measurements. This method consisted of simultaneously recovering stellar birth radii and the metallicity evolution in the Galactic disc with radius and time, |$\rm [Fe/H]$|(R,τ), by varying the metallicity slope at any given lookback time so that the resulting |$\rm \text{R}_\text{birth}$| distributions remained physically meaningful. It should be emphasized that this technique was almost completely independent of any analytical or numerical chemo-dynamical models; the only constraint used was that stars were required to be born at physically meaningful radii and the |$\rm \text{R}_\text{birth}$| distributions of older populations peaked at smaller radii, as expected for inside-out disc formation. Once in possession of |$\rm \text{R}_\text{birth}$|⁠, Minchev et al. (2018) constrained a number of chemo–kinematic relations with only a handful of stars from AMBRE:HARPS (De Pascale et al. 2014; Hayden et al. 2017) and HARPS-GTO (Adibekyan et al. 2012; Delgado Mena et al. 2017), combined with stellar ages computed with the starhorse code (Anders et al. 2018; Queiroz et al. 2018). Recently, Lu et al. (2022a) discovered the linear relation between metallicity scatter for a given age and the metallicity birth gradient evolution with cosmic time, thus making the |$\rm \text{R}_\text{birth}$| determination fully self-consistent. The results of the above work have opened up a new way of studying the chemo-dynamical evolution of the Milky Way, where instead of relying on forward modelling over cosmic time (cosmological simulation or semi-analytical modelling), one gets there directly from the data with few prior assumptions on the detailed chemical enrichment.

If |$\rm \text{R}_\text{birth}$| can be recovered using the observed age and |$\rm [Fe/H]$| (Minchev et al. 2018; Lu et al. 2022a), as it also appears to be the case in cosmological simulations of Milky Way-like galaxy discs (Feltzing, Bowers & Agertz 2020), then this should have implications for all other chemical elements. Indeed, Ness et al. (2019), Sharma et al. (2022), and Carrillo et al. (2023) (using APOGEE, galah, and cosmological zoom-in simulations, respectively) demonstrated that for stars in small bins of |$\rm [Fe/H]$| and age, all other abundances had very small scatter. More recently, Hayden et al. (2022) used the correlation between age, |$\rm [Fe/H]$|⁠, and elemental abundances to estimate ages for the galah dataset, thus using properties of stars born at the same time and place without actually recovering |$\rm \text{R}_\text{birth}$| (see also He, Luo & Chen 2022).

In this paper, we investigate the time evolution of chemical abundances across the Galactic disc using |$\rm \text{R}_\text{birth}$| derived in this work for APOGEE DR17 stars with ages derived by Anders et al. (2023) following the method of Lu et al. (2022a). Specifically, the questions we wish to answer are: (i) what do the α-sequences look like across time and |$\rm \text{R}_\text{birth}$|⁠, and (ii) how do abundance gradients evolve, and what do they reveal about the different nucleosynthetic processes? |$\rm \text{R}_\text{birth}$| estimates for a high resolution spectroscopic survey, like APOGEE, will provide light on these questions, as we want to probe the Galaxy’s evolution in as fine detail as possible. While |$\rm \text{R}_\text{birth}$| can only be inferred under some modelling assumptions, the relative |$\rm \text{R}_\text{birth}$| for a given age will be a useful tool in investigating the Milky Way’s past.

The paper is organized as follows. The data and method for estimating |$\rm \text{R}_\text{birth}$| used in this work are discussed in Sections 2 and 3, respectively. Section 4 details our results, with analyses of |$\rm [X/H]$| and |$\rm [X/Fe]$| abundance gradients across time in Sections 4.1 and 4.2, and the time evolution of the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$|⁠, age–|$\rm [X/Fe]$|⁠, and age–|$\rm [X/H]$| planes is discussed in Section 4.3. Sections 5 and 6 provide the discussion and key conclusions of this paper.

2 DATA

We use data provided by the 17th data release of APOGEE (Majewski et al. 2017; Abdurro’uf et al. 2022) of the fourth phase of the Sloan Digital Sky Survey (SDSS-IV; Blanton et al. 2017). The APOGEE survey collects data from both the northern and southern hemispheres using two high-resolution spectrographs (Wilson et al. 2019) on the 2.5 m Sloan Foundation telescope at Apache Point Observatory (Gunn et al. 2006) and the 2.5 m du Pont telescope at Las Campanas Observatory (Bowen & Vaughan 1973). The data reduction pipeline is outlined in Nidever et al. (2015), and stellar parameters and abundances are processed by the APOGEE Stellar Parameter and Chemical Abundance Pipeline (ASPCAP; Holtzman et al. 2015; García Pérez et al. 2016; Jönsson et al. 2020) and the APOGEE line list (Smith et al. 2021).

We partner the abundances with kinematics from the astroNN catalogue (Mackereth & Bovy 2018; Leung & Bovy 2019) and spectroscopic stellar ages for 185 282 stars from Anders et al. (2023), which have a median uncertainty of 1 Gyr and relative uncertainty <25 per cent. The ages are estimated from stellar parameters and elemental abundance ratios (most importantly T_eff, [C/Fe], [N/Fe], [Mg/Fe], and log g; not including [Fe/H]) using the supervised machine learning technique XGBoost trained on a set of more than 3000 red-giant and red-clump stars with asteroseismic ages (Miglio et al. 2021). We calculate the guiding radius (⁠|$\rm \text{R}_\text{guide}$|⁠) from the cylindrical rotational velocity:

$$\begin{eqnarray} {\rm \text{R}_\text{guide}}= {\rm R_\text{gal}}V_\phi /V_0, \end{eqnarray}$$

where V_ϕ is the Galactocentric azimuthal velocity and V₀ = 229.76 km s⁻¹ is the Milky Way rotation curve at solar radius (Schönrich, Binney & Dehnen 2010; Bovy et al. 2012b).

Since our |$\rm \text{R}_\text{birth}$| determination method is directly dependent on age and |$\rm [Fe/H]$|⁠, we select a sample with relatively small uncertainties (σ_age < 1.5 Gyr and |$\sigma _{\rm {\rm [Fe/H]}}\lt 0.015$| dex), and remove flagged |$\rm [Fe/H]$|⁠. In order to avoid systematic abundance trends with evolutionary state (Gruyters, Nordlander & Korn 2014; Dotter et al. 2017; Souto et al. 2018; Jofré, Heiter & Soubiran 2019; Liu et al. 2019; Shetrone et al. 2019; Lian et al. 2022), we focus on red giant (⁠|$2 \lt \rm{\log \,\,{\it g}}\lt 3.6$|⁠, and 4250 K |$\le {T_{\rm eff}}\le 5500$| K) disc (|z| ≤ 1 kpc, eccentricity <0.5, |$|{\rm [Fe/H]}| \lt 1$|⁠) stars. The left panel of Fig. 1 illustrates the log g and T_eff cuts employed in this work. This leaves us with 145 447 stars to use for inferring birth radii.

Figure 1.

Left: 2D histogram of the Kiel diagram of the APOGEE DR17 catalogue with Anders et al. (2023) ages. The black dashed lines indicate the region of the red giant branch that is used in this work, which is chosen to minimize the effects of atomic diffusion and first dredge-up on abundances and ages. Middle: |$\rm R_\text{gal}$|–z spatial distribution of APOGEE DR17 red giant stars with abundance and age quality cuts as described in Section 2. The black dashed lines show the z cuts we used to choose disc stars. Right: 2D histogram of the Galactic x–y plane of the 87 426 red giant disc stars used to measure abundance variation across time. The sun is located at (8.125, 0, 0.02) kpc (GRAVITY Collaboration 2018; Bennett & Bovy 2019). Only bins with >10 stars are shown.

For the abundance gradient analysis, we use calibrated ASPCAP APOGEE DR17 abundances [|$\rm [X/Fe]$| for X  = Mg, Ca (α), Mn (iron-peak), Al (light odd Z), C (light), Ce² (s-process)]. To ensure the highest quality sample, we only keep stars with unflagged abundances, and |$\rm [X/Fe]$| measurement uncertainties <0.1 dex for Ce and <0.03 dex for the other elements. We also remove spurious measurements using stars with |$-1 \le {\rm [X/Fe]}\le 1$|⁠, and employ a signal-to-noise ratio cut of SNREV > 65.

Due to the few stars with ages >12 Gyr, we are not able to reliably estimate the metallicity gradient beyond then.³ Therefore, we focus our analysis on stars with age ≤12 Gyr. This gives us a sample size of 87 426 stars to measure how abundances vary across the disc with cosmic time. The distribution of these stars in the Galactic x–y plane is shown in the right panel of Fig. 1.

3 METHOD FOR ESTIMATING BIRTH RADII

We follow the method from Lu et al. (2022a) to estimate birth radii using stellar ages and metallicities. This method assumes that the birth metallicity gradient is always linear in radius, as supported by both observations (Deharveng et al. 2000; Esteban et al. 2017; Arellano-Córdova et al. 2021) and simulations (Vincenzo & Kobayashi 2018; Lu et al. 2022b). Therefore for any lookback time (τ) we can write |$\rm [Fe/H]$| as a function of the metallicity gradient at that time (⁠|$\rm \nabla [Fe/H](\tau)$|⁠), birth radius, and the metallicity at the Galactic centre (⁠|$\rm [Fe/H](R = 0, \tau$|⁠)):

$$\begin{eqnarray} {\rm [Fe/H]}(\rm{\rm \text{R}_\text{birth}}, \tau) = \rm{\rm \nabla [Fe/H](\tau)}\rm{\rm \text{R}_\text{birth}}+ \rm{\rm [Fe/H](R = 0, \tau )}. \end{eqnarray}$$

(1)

Rearranging equation (1), we are able to estimate |$\rm \text{R}_\text{birth}$| as a function of metallicity and age:

$$\begin{eqnarray} \rm{\rm \text{R}_\text{birth}}(\text{age}, {\rm [Fe/H]}) = \frac{{\rm [Fe/H]}-\rm{\rm [Fe/H](R = 0, \tau )}}{\rm{\rm \nabla [Fe/H](\tau)}}, \end{eqnarray}$$

(2)

where τ is lookback time and age is observable. This method also assumes the star-forming gas in the Milky Way is azimuthally chemically homogeneous, in agreement with observations (Arellano-Córdova et al. 2021; Esteban et al. 2022), and that the full intrinsic metallicity range is well sampled in each age bin.

The time evolution of the metallicity at the Galactic centre, |$\rm [Fe/H](R = 0, \tau$|⁠), is estimated from the upper envelope of the sample’s age–metallicity relation (AMR) for ≥10 Gyr. The upper envelope is defined as the maximum |$\rm [Fe/H]$| for a given age bin, after removing outliers using the interquartile method. The smooth enrichment of the Galactic centre for τ < 10 Gyr is modelled by a log function, where |$\rm [Fe/H]$| (R = 0, τ = 0) is chosen under the assumptions that the youngest stars near the solar neighbourhood have a metallicity of 0.06 dex (Nieva & Przybilla 2012; Asplund, Amarsi & Grevesse 2021) and the present-day gradient is −0.074 dex kpc⁻¹ for the youngest stars (Zhang et al. 2021; Myers et al. 2022).

Using the cosmological hydrodynamic simulations of Buck et al. (2019); Libeskind et al. (2020); Khoperskov et al. (2022b), Lu et al. (2022a) showed that the range of [Fe/H] (defined as the 95th minus the fifth percentile) in a given age bin over time represents the shape of the metallicity gradient as a function of lookback time:

$$\begin{eqnarray} \rm{\rm \nabla [Fe/H](\tau)}= a\, \text{Range}\widetilde{{\rm [Fe/H]}}(age) + b, \end{eqnarray}$$

(3)

where parameters a and b are constants estimated from the boundary conditions and |$\text{Range}\widetilde{{\rm [Fe/H]}}(age)$| is the metallicity range normalized to have a minimum of 0 and a maximum of 1 across all age bins. When the scatter in |$\rm [Fe/H]$| is the smallest across all ages (⁠|$\text{Range}\widetilde{{\rm [Fe/H]}}(age) =0$|⁠), |$\rm{\rm \nabla [Fe/H](\tau)}= b$|⁠, which sets b as the current gradient of star-forming gas in the Milky Way, or −0.074 dex kpc⁻¹ (Zhang et al. 2021). When the scatter in [Fe/H] is the largest across all ages (⁠|$\text{Range}\widetilde{{\rm [Fe/H]}}(age) =1$|⁠), |$\rm{\rm \nabla [Fe/H](\tau)}= a + b$|⁠. This value is the steepest we would expect the metallicity gradient to have ever been, which we set to −0.15 dex kpc⁻¹ as found in other works (Minchev et al. 2018; Lu et al. 2022a), and thus forces a = −0.076 dex kpc⁻¹. Different values for the steepest metallicity gradient (and therefore a) are explored in Appendix  C. After removing outliers using the interquartile method, we employ the standard deviation in |$\rm [Fe/H]$| as a proxy for the metallicity range for 2 < age < 12 Gyr. The present-day standard deviation in |$\rm [Fe/H]$| is then extrapolated using a linear fit from the 2–6 Gyr values. Finally, to estimate |$\rm \text{R}_\text{birth}$| we plug equation (3) into the denominator of equation (2). The average uncertainty in |$\rm \text{R}_\text{birth}$| is 0.7 kpc, estimated from 100 Monte Carlo samples.

The left panel of Fig. 2 shows the measured scatter in metallicity (blue curve) and the inferred radial gradient (red curve). To best trace the substructure in the data, we estimate |$\rm [Fe/H](R = 0, \tau$|⁠) and |$\rm \nabla [Fe/H](\tau)$| using age bins of width 0.5 Gyr and smooth the functions over a 1-Gyr window (i.e. calculate the running mean with a window size of 3). This ensures a smooth metallicity gradient that removes artificial discontinuities and minimizes artefacts due to noise. We find a quick steepening in |$\rm \nabla [Fe/H](\tau)$| between 12 and 10 Gyr ago, which Lu et al. (2022a) linked to the effect of the Milky Way last massive merger – the Gaia–Sausage–Enceladus (GSE) event (Belokurov et al. 2018; Helmi et al. 2018). Instead of the monotonic flattening of the metallicity slope for age < 9 Gyr found by Lu et al. (2022a), however, with our dataset we recover two additional fluctuations in |$\rm \nabla [Fe/H](\tau)$| ∼6 and ∼4 Gyr ago.

Figure 2.

Left: The scatter in |$\rm [Fe/H]$| for a given age bin (blue) and the inferred metallicity gradient (red). The standard deviation in metallicity and the gradient are related through equation (3), once the |$\rm [Fe/H]$| standard deviation is normalized between 0 and 1. Right: The AMR of 145 447 APOGEE DR17 red-giant disc stars with ages from Anders et al. (2023). The colour-coded lines represent the running mean of mono-|$\rm \text{R}_\text{birth}$| populations with Δ|$\rm \text{R}_\text{birth}$|=0.75 kpc and mean values given by the colour bar. The shape of each |$\rm \text{R}_\text{birth}$| track is determined from |$\rm [Fe/H](R = 0, \tau$|⁠) and |$\rm \nabla [Fe/H](\tau)$|⁠. The fluctuations in |$\rm \nabla [Fe/H](\tau)$| at ∼6 and ∼4 Gyr seen in the left panel break the monotonic increase in |$\rm [Fe/H]$| for outer-disc mono-|$\rm \text{R}_\text{birth}$| populations and even appear to dilute the ISM 4–5 Gyr ago.

The right panel of Fig. 2 presents the AMR stellar density distribution of our APOGEE red giant sample. The colour-coded curves indicate the running mean of mono-|$\rm \text{R}_\text{birth}$| populations. The fast increase in metallicity between 12 and 10 Gyr ago corresponds to the quick steepening in the gradient, as seen in the left panel. This initial burst in star formation has long been linked to the formation of the high-α population in the |$\rm [Fe/H]$|–|$\rm [\alpha /Fe]$| plane (e.g. Chiappini et al. 1997; Matteucci 2021). The strong ripples in |$\rm \nabla [Fe/H](\tau)$| at lookback times of ∼6 and ∼4 Gyr seen in the left panel interrupt the monotonic increase in |$\rm [Fe/H]$| for mono-|$\rm \text{R}_\text{birth}$| populations at larger radii, and even appear to dilute the ISM in the range 4 ≲ age ≲ 5 Gyr. We discuss the possible causes of these fluctuations in Section 5.2.

Fig. 3 contrasts the difference between |$\rm \text{R}_\text{birth}$| (top) and |$\rm \text{R}_\text{guide}$| (bottom) for solar neighbourhood stars. Since |$\rm \text{R}_\text{birth}$| is not directly observable, we validate our estimates by comparing to our current understanding of how mono-age |$\rm \text{R}_\text{birth}$| distributions should behave. In agreement with expectations, we find that the youngest stars peak locally, as they have had less exposure to radial mixing processes (e.g. Minchev et al. 2013; Feltzing et al. 2020; Netopil et al. 2022), while older population peaks shift to progressively lower radii as expected for inside-out disc formation. Our method also reproduces |$\rm \text{R}_\text{birth}$| distributions found in state-of-the-art simulations (e.g. Agertz et al. 2021; Carrillo et al. 2023). In contrast to Minchev et al. (2018), who imposed such shifting in the peaks of mono-|$\rm \text{R}_\text{birth}$| populations in order to estimate the Milky Way’s |$\rm \nabla [Fe/H](\tau)$|⁠, here we reach a very similar configuration directly from the data, as described above.

Figure 3.

Top: |$\rm \text{R}_\text{birth}$| and bottom: |$\rm \text{R}_\text{guide}$| density distributions for mono-age populations in the solar neighbourhood (⁠|$7.9 \lt \rm{\rm R_\text{gal}}\lt 8.35$| kpc, |z| <.3 kpc), with each conditional density scaled by the total number of observations such that the total area under all the densities sums to 1. Within the solar neighbourhood, we find stars born throughout the disc that have migrated here over time. The youngest stars are found to be born locally, while the oldest stars formed further inwards. This is in agreement with expectations, as young stars have had less time to migrate, while older ones have had longer exposure to radial migration processes.

4 RESULTS

With the birth properties of stars, namely their ages and birth radii, we can record the Milky Way disc’s evolution with lookback time. Using |$\rm \text{R}_\text{birth}$| derived by the method in Section 3, we start our analysis by investigating the birth radial abundance gradients in |$\rm \text{R}_\text{birth}$|–|$\rm [X/H]$| (Section 4.1) and |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| (Section 4.2) planes to understand the overall enrichment of the Galactic disc. In Section 4.3, we then explore the time evolution of mono-|$\rm \text{R}_\text{birth}$| populations in the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$| and age–|$\rm [X/Fe]$| abundance planes to understand the relative enrichment of different nucleosynthetic sources.

4.1 [X/H] variation with |$\rm \text{R}_\text{birth}$|

We begin by looking at an array of elements to examine how the different nucleosynthetic families evolved throughout the Galactic disc. The variations of |$\rm [X/H]$| with |$\rm \text{R}_\text{birth}$| for mono-age populations, as colour-coded, are plotted in the top panels of Fig. 4. As expected, all |$\rm [X/H]$| trends are negative, similar to |$\rm [Fe/H]$|⁠. The effect of the |$\rm \nabla [Fe/H](\tau)$| fluctuation ∼9 Gyr ago can be readily seen in the oldest three age bins for all elements. The second and third fluctuations (∼4 and ∼6 Gyr) are not as apparent, although an overlap can be observed for the bins spanning the range 3–5 and 6–8 Gyr. For all elements and mono-age populations (except for Ce in the younger age bins), the relationships between |$\rm [X/H]$| and |$\rm \text{R}_\text{birth}$| could be fit relatively well with a linear model, and the slope variations are shown in the left panel of Fig. 5. Similar to predictions of both semi-analytical (e.g. Prantzos & Boissier 2000) and chemo-dynamical (e.g. Vincenzo & Kobayashi 2018) Milky Way-like models, we find that, overall, the slope of |$\rm [X/H]$| steepens with increasing lookback time. Most elements show multiple fluctuations in the |$\rm [X/H]$| gradient trends over lookback time (left panel of Fig. 5), similar to |$\rm \nabla [Fe/H](\tau)$|⁠. A notable exception is [Ce/Fe], whose gradient shows a continued steepening with increasing lookback time ∼6 Gyr ago, rather than inversion. This can be attributed to the strong overall flattening in its gradient with cosmic time (top rightmost panel of Fig. 4), which counteracts the two fluctuations ∼9 and 6 Gyr ago, seen well in the other elements.

Figure 4.

Time evolution of |$\rm [X/H]$| abundance gradients for our sample. Top: Running means of different mono-age populations overlaying the 87 426 red giant APOGEE DR17 stars in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/H]$| plane, determined by calculating the average |$\rm [X/H]$| for |$\rm \text{R}_\text{birth}$| bins of width 0.5 kpc, and then smoothed over 1 kpc. Bottom: Same as the top panels, but looking at the variation over |$\rm \text{R}_\text{guide}$|⁠. The total sample’s running mean is provided in each |$\rm \text{R}_\text{guide}$|–|$\rm [X/H]$| plane as a black dashed line. The |$\rm [X/H]$| gradient with lookback time has been significantly affected by radial migration, especially for the oldest mono-age populations. The younger populations tend to follow the non-conditioned running mean, indicating that there is significant scatter in |$\rm [X/H]$| about the mono-age running means. The |$\rm \text{R}_\text{birth}$|– and |$\rm \text{R}_\text{guide}$|–[Ce/H] trends using Ce derived in the BAWLAS catalogue are shown in Fig. B1 in the Appendix.

Figure 5.

Left: Measured |$\rm [X/H]$| radial gradients across lookback time, for Δlookback time  = 0.5 Gyr. The black dashed line is the derived |$\rm \nabla [Fe/H](\tau)$| from the left panel of Fig. 2. Most abundances show similar trends in |$\rm \text{R}_\text{birth}$|–|$\rm [X/H]$| with time, with slight variations caused by the rate and strength of different enrichment sources. Middle: Measured |$\rm [X/Fe]$| radial gradients over lookback time assuming a linear fit for each mono-age population. The maximum error on both |$\rm [X/H]$| and |$\rm [X/Fe]$| radial gradients is <0.005 dex kpc⁻¹, estimated from 100 Monte-Carlo samples where new abundances are drawn from normal distributions with the observed |$\rm [X/Fe]$| as the mean and their respective uncertainties as the standard deviation. Right: The measured |$\rm [X/Fe]$| radial gradient over lookback time (Δlookback time  = 1 Gyr) for |$\rm [Mn/Fe]$|⁠, |$\rm [C/Fe]$|⁠, and |$\rm [Ce/Fe]$| before (solid line) and after (dashed line) the respective break points (R_break; see Fig. 6). Mono-age populations show non-linear trends between these three abundances and |$\rm \text{R}_\text{birth}$|⁠, therefore we a fit a linear function before and after a break point, which is different for different ages and |$\rm [X/Fe]$| (Fig. 6). The radial |$\rm [Mn/Fe]$| and |$\rm [Ce/Fe]$| gradients are steeper before the break points, while |$\rm [C/Fe]$| sees an inversion in slope before and after the break point.

We contrast the |$\rm [X/H]$| variation with |$\rm \text{R}_\text{birth}$| to the variation with |$\rm \text{R}_\text{guide}$| in the bottom row of Fig. 4. Drastic differences can be seen, especially for the oldest mono-age bins, as also expected from chemo-dynamical simulations. The oldest populations show significantly flatter relations with |$\rm \text{R}_\text{guide}$|⁠, which is expected since older populations have been exposed longer to radial migration processes. The younger populations are most similar to the recovered |$\rm [X/H]$| radial birth gradients, though with larger scatter. The dashed black curve in each bottom panel shows the running mean of the total stellar sample, whose slope and shape is strongly dependent on the survey selection function and the quality cuts used for the sample at hand and in general does not provide any useful information. Most of the mono-age populations (age ≲8 Gyr) follow the same |$\rm \text{R}_\text{guide}$|–|$\rm [X/H]$| trend as the running mean of our sample while the mono-age populations separate in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/H]$| plane, indicating that there is minimal |$\rm \text{R}_\text{birth}$| information contained in |$\rm \text{R}_\text{guide}$|⁠.

4.2 [X/Fe] variation with |$\rm \text{R}_\text{birth}$|

We now explore the |$\rm [X/Fe]$| trends with |$\rm \text{R}_\text{birth}$|⁠, resulting from our recovered |$\rm [Fe/H]$| gradient evolution. The top panels of Fig. 6 show the evolution of the relative enrichment between Supernovae Ia contribution (as traced by Fe) and other sources over time throughout the Milky Way disc. The middle panel of Fig. 5 shows the value of each |$\rm [X/Fe]$| gradient as a function of lookback time.

Figure 6.

Same as Fig. 4, but for the |$\rm [X/Fe]$| abundance gradients. Top: Running means of mono-age populations in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| plane. Remarkably, in all elements except for Ce, a bimodality akin to that found in the |$\rm [Fe/H]$|–|$\rm [\alpha /Fe]$| plane can be seen, with the high-α sequence appearing as a distinct elongated blob composed of the oldest 2–3 age bins. This time period corresponds exactly to the period during which the initial steepening in |$\rm \nabla [Fe/H](\tau)$| takes place in Fig. 2. Most radial mono-age |$\rm [X/Fe]$| gradients can be adequately fit by a linear model, except for |$\rm [C/Fe]$|⁠, |$\rm [Ce/Fe]$|⁠, and |$\rm [Mn/Fe]$|⁠. In each mono-age population for these elements, we mark the location where the slopes either flatten or invert by a black dot with a white outline. Bottom: Same as top but using the guiding radius. Trends with both radius and time seen in the top panels are almost completely erased due the angular momentum redistribution in the disc that has taken place. The effect is strongest for the oldest populations, resulting in the trends with radius being completely lost and inverted for Mg, Ca, Mg, and Al. See Appendix Fig. B1 for the |$\rm \text{R}_\text{birth}$|– and |$\rm \text{R}_\text{guide}$|–|$\rm [Ce/Fe]$| relationship using |$\rm [Ce/Fe]$| from the BAWLAS catalogue.

Each age bin in Fig. 6 has its own unique trend between |$\rm [X/Fe]$| and |$\rm \text{R}_\text{birth}$| that does not necessarily capture the overall trend of the whole data. Specifically, starting with |$\rm [C/Fe]$| (top left panel of Fig. 6), we see that the full sample has a negative correlation between |$\rm [C/Fe]$| and |$\rm \text{R}_\text{birth}$|⁠. However, once conditioned on age, each mono-age population has a primarily weak |$\rm [C/Fe]$| radial birth gradient. The global near-zero radial gradient in |$\rm [C/Fe]$| over lookback time is also seen in the middle panel of Fig. 5, and may be a consequence of carbon’s relationship with mass and age (Salaris et al. 2015), and the relationship between Fe, age, and birth radius. Since |$\rm [C/Fe]$| has a minimal radial gradient with age, the [C/H] gradient varies most similar to |$\rm [Fe/H]$| over time (left panel of Fig. 5).

The radial |$\rm [Mg/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠, and |$\rm [Al/Fe]$| birth gradients for (almost) all mono-age populations are positive and flattening with time as seen in Fig. 6. Ca and Mg, being α-elements, show very similar trends in both their |$\rm [X/Fe]$| and [X/H] radial gradient evolution with time, except for the smaller range in |$\rm [Ca/Fe]$|⁠. Mg and Al also have essentially the same radial gradients (Fig. 5).

Unlike the α-elements, the |$\rm [Mn/Fe]$| radial gradient is negative, and has stayed at a fairly constant rate throughout time. Despite the fact that Mn, Mg, Ca, and Al have different enrichment sources, the 3–4 and 4–5 Gyr age tracks and 6–7 and 7–8 Gyr age tracks for these four elements live in very similar areas in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| plane, reflecting the |$\rm \nabla [Fe/H](\tau)$| fluctuations at those times (Fig. 6).

Another element ratio that stands out is |$\rm [Ce/Fe]$| showing a strong positive correlation with |$\rm \text{R}_\text{birth}$| for age <8 Gyr and a weaker relationship in older aged populations, including an inversion for oldest age bin. Unlike the other abundances, the age tracks do not cleanly separate in the |$\rm \text{R}_\text{birth}$|–|$\rm [Ce/Fe]$| plane. The weak |$\rm [Ce/Fe]$| abundance gradient for older aged stars suggests that the production rate of Ce and Fe varied slowly at first, then saw larger differences at later times (consistent with fig. 7 of Battistini & Bensby 2016).

So far, we have focused on the |$\rm [X/Fe]$| abundance ratios’ radial gradients in Fig. 6 being linear for each mono-age bin (Fig. 5, middle). However, some of the trends with |$\rm \text{R}_\text{birth}$| cannot accurately be fit with a linear model. Most noticeably, |$\rm [C/Fe]$|⁠, |$\rm [Ce/Fe]$|⁠, and the younger aged populations in |$\rm [Mn/Fe]$| have a break point, where the abundance trends begin to flatten, or even invert slope. The locations of the break points are shown in Fig. 6 as black dots with a white outline. The right panel of Fig. 5 shows these three abundances fit linearly before (solid line) and after (dashed line) their respective break points. The younger mono-age populations show an inverse in slope before and after the break point in |$\rm [C/Fe]$| and |$\rm [Ce/Fe]$|⁠, with a positive and negative radial gradient in the outer region of the disc, respectively. Meanwhile, the |$\rm [Mn/Fe]$| radial gradient is negative, with the gradient after the break point being flatter. These dramatic changes in the |$\rm \text{R}_\text{birth}$| gradients of |$\rm [C/Fe]$|⁠, |$\rm [Ce/Fe]$|⁠, and |$\rm [Mn/Fe]$| for different radii ranges should be considered in modelling the evolution of these elements.

We also show the trends of |$\rm [X/Fe]$| with guiding radius for different mono-age populations in the bottom panels of Fig. 6 to illuminate the additional information learned with knowledge of |$\rm \text{R}_\text{birth}$|⁠. Overall, we can see that many birth gradients are lost or altered as a result of radial migration. For instance, most |$\rm \text{R}_\text{guide}$|–|$\rm [X/Fe]$| relations are flatter or inverted, especially for the older age bins, compared to radial gradients measured with |$\rm \text{R}_\text{birth}$| (e.g. |$\rm [Al/Fe]$|⁠, |$\rm [Mg/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠), and non-linear gradients can become smoothed by stellar migration (e.g. |$\rm [C/Fe]$|⁠).

Radial gradients with |$\rm \text{R}_\text{birth}$| also reveal structure in the data that recreate the chemical evolution of the Galaxy which cannot be seen in the present-day disc. The top panels of Fig. 6 show a bimodality in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| plane, particularly for the α-elements Mg and Ca, the odd-Z element Al, and the iron-peak element Mn. One of the groups is captured by older (>9 Gyr) stars, has high |$\rm [Mg/Fe]$|⁠, and is representative of the high-α sequence, while the other corresponds to the low-α sequence. The clean separation of the high- and low-α sequences is not typically seen in the |$\rm \text{R}_\text{guide}$|–|$\rm [X/Fe]$| plane (see bottom panels of Fig. 6). With a clear view of the different enrichment histories, we see that the high-α sequence typically has a steeper |$\rm [X/Fe]$| radial birth gradient, suggesting that it formed in a chemically non-homogeneous environment. Within the low-α sequence, the younger stars appear to separate into their own distinct group in the light element abundances and |$\rm [Mg/Fe]$|⁠. This is due to a rapid weakening of the |$\rm [X/Fe]$| gradient 3–4 Gyr ago.

4.3 Stellar chemical enrichment with lookback time

Until here, we have explored how the radial gradients of individual chemical abundances and abundance ratios evolve across the Galactic disc with time. In this section, we now study the time evolution of mono-|$\rm \text{R}_\text{birth}$| populations in the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$| and age–|$\rm [X/Fe]$| planes to investigate the relative rate of enrichment. As a caveat, we note that our age estimates, while being independent of [Fe/H], do depend to some degree on [Mg/Fe] (see Anders et al. 2023). However, this dependency seems to have a minimal effect on our results as we are able to reproduce the same trends using ages derived using isochrone fitting (see Section  A in the Appendix).

4.3.1 Time evolution of the [Fe/H]–[Mg/Fe] plane

We begin this section by looking at the formation and evolution of the α-sequences across |$\rm \text{R}_\text{birth}$| and age. Fig. 7 shows different |$\rm \text{R}_\text{birth}$| tracks as a function of time in the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$| plane, similar to Buck (2020), Sharma et al. (2021), and Lu et al. (2022a). For a given |$\rm \text{R}_\text{birth}$|⁠, the track begins with a higher value of |$\rm [Mg/Fe]$|⁠, and over time moves to lower values of |$\rm [Mg/Fe]$| and higher values of |$\rm [Fe/H]$| for younger aged stars. Overall, the evolution through this plane for a given |$\rm \text{R}_\text{birth}$| is smoothly increasing in |$\rm [Fe/H]$|⁠. There is, however, a wiggle at a lookback time of ∼3.5–5 Gyr where younger stars jump to lower |$\rm [Fe/H]$| at higher |$\rm \text{R}_\text{birth}$|⁠, and see minimal changes in |$\rm [Fe/H]$| at smaller |$\rm \text{R}_\text{birth}$|⁠. This non-monotonic evolution in the α-bimodality is caused by the fluctuation in |$\rm \nabla [Fe/H](\tau)$| at ∼4 Gyr, seen in Fig. 5, and leads to three different ages having almost the same |$\rm [Fe/H]$| for a given |$\rm \text{R}_\text{birth}$| (also seen in the right panel of Fig. 2). During this time of the |$\rm \nabla [Fe/H](\tau)$| fluctuation ∼6 Gyr ago, the |$\rm [Fe/H]$| and |$\rm [Mg/Fe]$| gradients stay fairly constant with time, causing stars born 5.5–8 Gyr ago to have similar |$\rm [Fe/H]$| and |$\rm [Mg/Fe]$| at a given |$\rm \text{R}_\text{birth}$|⁠. Note that the points along each |$\rm \text{R}_\text{birth}$| track are all separated by 1 Gyr, i.e. a smaller separation corresponds to a slower enrichment. This is another reason why mono-abundance populations cannot provide a good proxy for mono-age populations, in addition to radial mixing or sample covering larger Galactic radius (Minchev et al. 2017).

Figure 7.

The density distribution of 87 426 red giant APOGEE stars in the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$| plane with the time evolution of |$\rm \text{R}_\text{birth}$| tracks overlaid. The tracks are chosen with Δ|$\rm \text{R}_\text{birth}$| = 1 kpc and Δage  = 0.5 Gyr, with points at every other age bin. The fluctuation in the |$\rm [Fe/H]$| gradient at ∼9 Gyr corresponds to the valley between the high- and low-α sequences, while the fluctuation at ∼6 Gyr corresponds to a period of slower enrichment.

Fig. 7 reveals that the high-α sequence (defined in the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$| plane such as in Fig. F1) has formed close to the Galactic centre (⁠|$\rm{\rm \text{R}_\text{birth}}\le 7$| kpc), in agreement with Buck (2020) and Lu et al. (2022a) and as expected for an inside-out disc formation. This is also consistent with high-α stars being predominantly found in the inner parts of the Galaxy (e.g. Bensby et al. 2011; Anders et al. 2014; Hayden et al. 2015). It is uniformly old (yellow dots) up to |$\rm [Fe/H]$|≈0–0.1 dex and its age decreases to about 6.5 Gyr at the most metal-rich part (where the low- and high-α sequences meet). The majority of the stars composing the high-α sequence were formed at the very inner disc, according to the colour-coded curves. As found also by Lu et al. (2022a), the high-α sequence formed during the early-on steepening of |$\rm \nabla [Fe/H](\tau)$|⁠. Going through the steepest point at lookback time ∼9 Gyr caused a large jump in both |$\rm [Fe/H]$| and |$\rm [Mg/Fe]$|⁠, thus forming the gap in the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$|plane.

Fig. 8 shows the distribution of the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$| plane as a function of |$\rm \text{R}_\text{birth}$| (top) and |$\rm \text{R}_\text{guide}$| (bottom), allowing us to make direct comparisons between the observable and inherent trends. While the high-α sequence is present beyond the solar neighbourhood in |$\rm \text{R}_\text{guide}$|⁠, with a small proportion of stars currently located in the outer disc, it was predominantly formed at lower |$\rm \text{R}_\text{birth}$|⁠, with only about 2  per cent of its stars born beyond 7 kpc. The formation structure of the high-α sequence, which is typically concluded to be radially uniform in the inner disc and solar neighbourhood (Bensby et al. 2012; Nidever et al. 2014; Hayden et al. 2015), is also masked by radial migration. The most metal-rich high-α stars predominantly formed closer to the Galactic centre, while the lower |$\rm [Fe/H]$| stars were formed further out, suggesting an inside-out formation also for the high-α sequence.

Figure 8.

Distribution of stars in the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$| plane separated as a function of Top: birth radius and Bottom: guiding radius. The grey points are the individual stars within each |$\rm \text{R}_\text{birth}$|/|$\rm \text{R}_\text{guide}$| bin, with the filled in density contours on top. The spatial separation as a function of |$\rm \text{R}_\text{birth}$| creates tighter relationships between |$\rm [Mg/Fe]$| and |$\rm [Fe/H]$|⁠, especially in the low-α sequence. The white points and labels in the top row denote the lookback time in Gyr.

Unlike the high-α sequence, the low-α sequence (see Fig. F1 for our visual definition of the two sequences) was formed throughout the Galactic disc. Fig. 8 shows that the overall trend of the low-α sequence evolves with |$\rm \text{R}_\text{birth}$| similar to |$\rm \text{R}_\text{guide}$|⁠. That is, the low-α stars have lower values of |$\rm [Fe/H]$| for higher |$\rm \text{R}_\text{birth}$| or |$\rm \text{R}_\text{guide}$|⁠, indicating that the metal-rich low-α stars formed at very low |$\rm \text{R}_\text{birth}$| (consistent with suggestions by e.g. Kordopatis et al. 2015). Unlike with |$\rm \text{R}_\text{guide}$|⁠, separating by |$\rm \text{R}_\text{birth}$| creates much narrower distributions, revealing a tight relationship between |$\rm [Mg/Fe]$| and |$\rm [Fe/H]$| in the low-α sequence for a given |$\rm \text{R}_\text{birth}$|⁠. Specifically, we can see regions of higher density in this plane that show the transition of the location of the low-α sequence with Galactic birth location. This indicates that the higher |$\rm [Mg/Fe]$| tail of the low-α sequence formed first in the inner disc, while the stars composing the lower |$\rm [Mg/Fe]$|⁠, higher |$\rm [Fe/H]$| end formed later on further out in the disc. The unveiling of this spatial structure with |$\rm \text{R}_\text{birth}$| illustrates that the high-|$\rm [Fe/H]$| region of the low-α sequence migrated outwards, blurring the strong inherent correlations between |$\rm [Mg/Fe]$|⁠, |$\rm [Fe/H]$|⁠, and age.

4.3.2 Age–abundance plane

The age–|$\rm [X/Fe]$| plane provides insight into the overall evolution of |$\rm [X/Fe]$| with time. The |$\rm [X/Fe]$| abundance ratios of Mg, Ca, Mn, and Al have been fairly constant over the past ∼9 Gyr, with a distinct steeper slope at 9–12 Gyr (grey points in top panels of Fig. 9). Conditioning on mono-|$\rm \text{R}_\text{birth}$| populations changes this perception, as shown by the colour-coded curves. Except for |$\rm [C/Fe]$| – whose mono-|$\rm \text{R}_\text{birth}$| populations show minimal differences due to the already very strong positive age–|$\rm [C/Fe]$| relation – the mono-|$\rm \text{R}_\text{birth}$| populations occupy different |$\rm [X/Fe]$| for a given age, and show a strong age–|$\rm [X/Fe]$| correlation even when the overall population does not. Each mono-|$\rm \text{R}_\text{birth}$| population for |$\rm{\rm \text{R}_\text{birth}}\lt 9$| kpc increases in |$\rm [Mg/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠, and |$\rm [Al/Fe]$| until 4 Gyr, when the |$\rm [X/Fe]$| abundance increases at a slightly faster rate for ∼5 Gyr, before proceeding to enrich at a faster rate until 12 Gyr. An exception to this is that the innermost |$\rm \text{R}_\text{birth}$| see a decrease in |$\rm [Al/Fe]$| until about 8 Gyr ago, then proceed to increase in |$\rm [Al/Fe]$|⁠. The |$\rm{\rm \text{R}_\text{birth}}\gt 9$| kpc region is more concentrated towards younger ages, and therefore only sees the steep enrichment at age <6 Gyr. |$\rm [Mn/Fe]$| follows a similar pattern, just with a negative correlation with age. Unlike the other abundances, the correlation between |$\rm [Ce/Fe]$| and age is strongly dependent on the birth radius. For outer radii, |$\rm [Ce/Fe]$| declines with increasing age, while in the inner disc the relationship inverts. The above results illustrate the importance of incorporating |$\rm \text{R}_\text{birth}$| estimates in studying the Galaxy’s evolution, as in their absence it can be concluded that |$\rm [\alpha /Fe]$| is not a good tracer for age up to the slope increase at ∼9 Gyr, as recently concluded by Hayden et al. (2022).

Figure 9.

Mono-|$\rm \text{R}_\text{birth}$| populations in the top: age–|$\rm [X/Fe]$| and bottom: age–|$\rm [X/H]$| planes of 87 426 red giant APOGEE DR17 stars. The running means are derived using age bins of 0.5 Gyr and averaged over 1 Gyr. The running mean of the total population is given in the bottom panels, calculated the same way as the mono-|$\rm \text{R}_\text{birth}$| bins. The mono-|$\rm \text{R}_\text{birth}$| populations lie on steeper tracks in the age–|$\rm [X/Fe]$| and –|$\rm [X/H]$| planes than the overall population. Similar to |$\rm [Fe/H]$|⁠, |$\rm [X/H]$| has three fluctuations in the gradient, so a given |$\rm \text{R}_\text{birth}$| is not monotonically increasing in |$\rm [X/H]$| with lookback time. Note that when working in terms of mono-|$\rm \text{R}_\text{birth}$| populations, age becomes lookback time, i.e. corresponds to the actual evolution with cosmic time. See Fig. B1 in the Appendix for age–|$\rm [Ce/Fe]$| and –[Ce/H] trends using Ce from the BAWLAS catalogue.

The bottom panels of Fig. 9 show the age–|$\rm [X/H]$| plane, which is equivalent to the AMR shown in the right panel of Fig. 2. Similar to the age–|$\rm [X/Fe]$| plane, separating the age–|$\rm [X/H]$| plane by |$\rm \text{R}_\text{birth}$| provides structure to an otherwise uncorrelated relationship. The black dashed line represents the full sample’s running mean in |$\rm [X/H]$| across age. |$\rm [X/H]$| overall shows lots of scatter and minimal correlation with age, particularly for age <5 Gyr. Conditioning the sample on |$\rm \text{R}_\text{birth}$| creates structure similar to the age–|$\rm [Fe/H]$| plane, given in the right panel of Fig. 2. Each mono-|$\rm \text{R}_\text{birth}$| population shows a clear enrichment in |$\rm [X/H]$| with time, with the same fluctuations at 4, 6, and 9 Gyr (except perhaps for Ce).

With |$\rm \text{R}_\text{birth}$| available, it is clear from Fig. 9 that both |$\rm [X/H]$| and |$\rm [X/Fe]$| correlate well with age for all birth radii, were it not for radial migration blurring the relations. The statistical phenomenon that causes this flattening and even inversion of relations seen in Figs 6 and 9 is known as Yule–Simpson’s paradox (Simpson 1951), which was recently shown to be ubiquitous in the field of Galactic Archaeology (Minchev et al. 2019).

5 DISCUSSION

5.1 Comparison to radial gradients found in other work

Due to the angular momentum redistribution of stars within the Milky Way disc throughout its lifetime, the radial abundance gradients measured today for mono-age populations do not reflect those at birth. Observed present-day [X/H] (typically |$\rm [Fe/H]$| or [M/H]) and |$\rm [X/Fe]$| gradients are found to be flatter for older ages or show only weak evolution with age (Casagrande et al. 2011; Anders et al. 2017, 2023; Myers et al. 2022; Sales-Silva et al. 2022; Ratcliffe & Ness 2023). Using our recovered |$\rm \text{R}_\text{birth}$|⁠, we find a distinctly different gradient evolution with age, which in this case can only be thought of as lookback time.

The inherent birth gradients in |$\rm [X/Fe]$| that we recover are much steeper than the observed gradients of mono-age populations (>0.02 dex kpc⁻¹ for |$\rm [Ce/Fe]$| and |$\rm [Al/Fe]$|⁠, <−0.02 dex kpc⁻¹ for |$\rm [Mn/Fe]$|⁠, >0.02 dex kpc⁻¹ for lookback time >9 Gyr in |$\rm [Mg/Fe]$|⁠). The |$\rm [X/H]$| gradients (left panel of Fig. 5) show a clear flattening with time, in agreement with Prantzos & Boissier (2000); Vincenzo & Kobayashi (2018); Agertz et al. (2021); and Magrini et al. (2023). This is contradictory to some chemical evolution models and simulations, which assume the abundance gradients are steepening with cosmic time (e.g. Chiappini, Matteucci & Romano 2001; Vincenzo & Kobayashi 2020; Bellardini et al. 2022), and illustrates the importance of utilizing the |$\rm \text{R}_\text{birth}$| information now available.

Fig. 6 shows the radial |$\rm [X/Fe]$| information lost to radial mixing processes, and explains the mildly positive/negative |$\rm [X/Fe]$| gradients observed in the present-day disc (Carrera & Pancino 2011; Yong, Carney & Friel 2012; Anders et al. 2017; Myers et al. 2022). The similarities we find between the evolution of |$\rm [Al/Fe]$| and |$\rm [Mg/Fe]$| with |$\rm \text{R}_\text{birth}$| agrees with previous findings that Al behaves like an α-element in the age–|$\rm [X/Fe]$| plane (Spina et al. 2021) and highlights the relative production of Al in SNII is similar to Mg.

While we approximated the radial birth |$\rm [X/Fe]$| gradients as linear, they can best be explained as a combination of two pieces. We fit a linear model both before and after a break point in |$\rm [C/Fe]$|⁠, |$\rm [Ce/Fe]$|⁠, and |$\rm [Mn/Fe]$|⁠, which are the abundances we deem to be the most affected by a break (right panel of Fig. 5). For each abundance, we find that the location of the break in the |$\rm [X/Fe]$| gradient is mostly moving to larger radii over time, starting at ∼4 kpc 10+ Gyr ago and located >8 kpc for stars 0–1 Gyr. For the α-elements and Al, their |$\rm [X/Fe]$| radial gradients steepen slightly after the break. The radial gradient in |$\rm [Mn/Fe]$| weakens after the break, and the break creates an inversion in slope for |$\rm [C/Fe]$| and |$\rm [Ce/Fe]$|⁠. The location of the break point varies for different |$\rm [X/Fe]$|⁠, and should be considered in modelling the chemical evolution of the Milky Way (e.g. a radially varying initial mass function; Guszejnov, Hopkins & Graus 2019; Horta et al. 2022).

5.2 What causes the fluctuations in |$\rm \nabla [Fe/H](\tau)$|?

As discussed in Section 3, we find three fluctuations in the metallicity gradient at ∼4, ∼6, and ∼9 Gyr ago. Although more work is needed to understand this, the older one may be attributed to the effect of the last massive merger, GSE (Belokurov et al. 2018; Helmi et al. 2018), as also concluded by Lu et al. (2022a) and supported by a very similar |$\rm \nabla [Fe/H](\tau)$| fluctuation seen in a Milky Way-like cosmological simulation by Buck et al. (2019) during a GSE-like merger (fig. A1 in Lu et al. 2022a). In agreement with this picture, Ciucă et al. (2022) found a dip and a stand-alone feature (the ‘blob’) in the AMR at the same relative time we find the first steepening in |$\rm \nabla [Fe/H](\tau)$|⁠, which they also attributed to the GSE. The recent study of Vincenzo et al. (2019) suggested that the GSE merger heated up the gas into the Milky Way halo, thus delaying star formation and forming the gap in the |$\rm [Fe/H]$|–|$\rm [\alpha /Fe]$| plane. The pause in star formation due to the merger would primarily be felt in the outer disc (seen in simulations; Buck et al. 2023), therefore the metallicity gradient steepens due to the continually enriching inner disc. This can be seen in Fig. 4, where the three oldest mono-age populations have similar |$\rm [X/H]$| at higher |$\rm \text{R}_\text{birth}$|⁠, but show an enrichment over time at the lowest |$\rm \text{R}_\text{birth}$|⁠.

We find the fluctuations in |$\rm \nabla [Fe/H](\tau)$| to be robust to both the range and standard deviation of different age bin widths. When accounting for age and metallicity uncertainty by sampling from a normal distribution with a mean of the observed value and a standard deviation of the measurement error, we find that the fluctuations at ∼4 and 6 Gyr ago become weaker. However, the non-monotonic flattening of the metallicity gradient is still seen (Fig. 10). We are also able to recover multiple fluctuations in the AMR using the age catalogues of Queiroz et al. (2023) (APOGEE DR17, GALAH DR3), Mackereth et al. (2019) (APOGEE DR17), Leung et al. (2023) (APOGEE DR17), and Sharma et al. (2018) (GALAH DR3); see also Lin et al. (2020) (GALAH DR2) and fig. 2 in Sahlholdt, Feltzing & Feuillet (2022). While Lu et al. (2022a) did not find multiple fluctuations in their derived |$\rm \nabla [Fe/H](\tau)$| using LAMOST, we find that performing a stricter cut in [Fe/H] error on their sample recovers a second fluctuation, though the fluctuations are weaker than recovered in the other surveys/data sets we examined. Investigating this further by performing a crossmatch between LAMOST and GALAH, we determine that LAMOST’s lower resolution may be the reason any additional structure beyond the single fluctuation at the time of GSE is concealed.

Figure 10.

The normalized scatter in |$\rm [Fe/H]$| (blue) and inferred |$\rm \nabla [Fe/H](\tau)$| (red) averaged over 100 Monte Carlo samples where both the age and |$\rm [Fe/H]$| are drawn from a normal distribution about their measured values and a standard deviation of their measurement uncertainty. To show the robustness of the fluctuations but still ensure a smooth evolution over time, we choose an age bin width of 1 Gyr over a rolling window of 0.2 Gyr (i.e. we calculate the scatter in |$\rm [Fe/H]$| for 1.5–2.5 Gyr, 1.7–2.7 Gyr, etc.).

It is tempting to attribute the |$\rm \nabla [Fe/H](\tau)$| fluctuations at ∼4 and 6 Gyr to passages by the Sagittarius dwarf galaxy (Ibata, Gilmore & Irwin 1994; Law & Majewski 2010; Laporte et al. 2018), at which time this minor (yet quite massive) merger is expected to have deposited |$\sim 50~{{\ \rm per\ cent}}$| of its gas to its host, thus diluting the disc in |$\rm [Fe/H]$| (Tepper-García & Bland-Hawthorn 2018). The effect of the gas infall from Sagittarius is expected to be the most prominent at the outskirts, diluting the gas metallicity in the disc, and thus directly affecting the radial metallicity gradient of stellar populations formed soon after close passages (Annem & Khoperskov 2022). Indeed, for both of these fluctuations, the metallicity gradient begins to steepen at the same time as Sagittarius’ estimated pericentre passages in literature (Law & Majewski 2010). Mor et al. (2019) found a period of star formation rate enhancement ∼2–3 Gyr ago, and suggested that it was possibly due to a recent merger. Indeed, we find a fluctuation just prior to this at ∼4 Gyr ago, which suggests that the star formation rate enhancement began after the effects of the merger on the metallicity gradient. The fluctuation ∼6 Gyr ago also coincides with the oldest of the three episodes of enhanced star formation (5.7 Gyr ago) found by Ruiz-Lara et al. (2020), which, in agreement with some recent simulations (Di Cintio et al. 2021; Khoperskov et al. 2022a), can be associated with a pericentre passage of Sagittarius.

While the |$\rm \nabla [Fe/H](\tau)$| wiggle at ∼9 Gyr does not show a strong signature in the AMR of mono-|$\rm \text{R}_\text{birth}$| populations (Fig. 2, right), a stagnation or decline with age is found in |$\rm [Fe/H]$|⁠, as well as in all other elements besides [Ce/H] where the effect is only weak, outside |$\rm \text{R}_\text{birth}$|∼6 kpc during the second and third events (Fig. 9, bottom). These stronger signatures in the AMR are due to the slower |$\rm [Fe/H]$| evolution in the innermost disc, after the initial star-formation burst ending around 8 Gyr ago in our data.

Using chemical abundance data from APOGEE and Gaia RVS (Gaia Collaboration 2022), Spitoni et al. (2023) suggest a three-infall chemical evolution model (albeit not taking into account radial migration), with the most recent gas infall starting ∼2.7 Gyr ago. While our three fluctuations in |$\rm \nabla [Fe/H](\tau)$| are synonymous with the gas infalls found by Spitoni et al. (2023), the most recent one affecting the volume probed by the APOGEE DR17 data was ∼4 instead of ∼2.7 Gyr ago.

5.3 The high-α sequence formed throughout the inner disc in a non-mixed environment

We find that the high-α sequence primarily formed within |$\rm{\rm \text{R}_\text{birth}}\le 7$| kpc (in agreement with state-of-the-art simulations, e.g. Buck 2020), with few high-α stars born at larger radii. This is resonant with Bensby et al. (2011); Anders et al. (2014); Hayden et al. (2015), and Queiroz et al. (2020), who find few high-α stars located further beyond the solar neighbourhood in the present-day Galactic disc, and it suggests that these stars have not had time to migrate further. The few high-α stars we find born at larger |$\rm \text{R}_\text{birth}$| are younger than the remainder of the high-α sequence, and in fact we find a mean difference in |$\rm \text{R}_\text{birth}$| of 3 kpc between the high-α sequence as a whole and the young (<4 Gyr) α-enhanced stars. The existence of these apparently young α-enhanced stars (Chiappini et al. 2015; Martig et al. 2015) has convincingly been suggested to be an artefact of incorrect age assignment due to binary mass transfer (Jofré et al. 2016, 2023; Fuhrmann & Chini 2017; Izzard et al. 2018). This suggests that their |$\rm \text{R}_\text{birth}$| estimates are incorrectly assigned, and explains the young α-enhanced stars’ inconsistent |$\rm \text{R}_\text{birth}$| with the remainder of the high-α sequence. Using the chemical ages of Anders et al. (2023), we find a ∼3 per cent contamination rate of young (<6 Gyr) stars in the high-α sequence, and therefore conclude they do not affect the results found in this work.

With knowledge of |$\rm \text{R}_\text{birth}$|⁠, we find that the high-α sequence is both a sequence of |$\rm \text{R}_\text{birth}$| and age, where the high-|$\rm [Fe/H]$| end is formed closer to the Galactic centre, and the low-|$\rm [Fe/H]$| end is formed further out (Figs 7 and 8). The high-α sequence began forming along steep |$\rm [Fe/H]$| and |$\rm [Mg/Fe]$| radial gradients, and does not follow one evolutionary track as concluded in some chemical evolution models and works (e.g. Haywood et al. 2015, 2019; Sharma et al. 2021). Given that (⁠|$\rm [Fe/H]$|⁠, |$\rm [Mg/Fe]$|⁠) correspond to different |$\rm \text{R}_\text{birth}$| tracks, this suggests that the gas was not chemically homogeneous during the formation of the high-α sequence (contrasting with what is observed today e.g. Eilers et al. 2022; Ratcliffe & Ness 2023), and in fact we find that the high-α sequence had larger |$\rm [X/Fe]$| gradients than the low-α sequence (in agreement with Buck et al. 2023). While previous works concluded that the high-α sequence evolved along a single chemical track due to efficient gas-phase mixing (e.g. Haywood et al. 2013, 2015; Nidever et al. 2014; Katz et al. 2021), the differences we see between the |$\rm \text{R}_\text{birth}$| tracks highlight the information lost to radial mixing processes (see Fig. 6).

5.4 The evolution of a given birth radius in the [Fe/H]–[Mg/Fe] plane with time

Fig. 7 reveals that the high-α sequence formed >10 Gyr ago⁴, which is the period during which the gradient is found to steepen (Fig. 2, left). The valley between the two sequences formed between 9 and 10 Gyr ago with the gap larger for lower |$\rm [Fe/H]$|⁠, captured by the steep gradients in [X/H]. The origin of the transition between the high- and low-α sequences is a fundamental question still not well understood. Recent work has proposed a variety of situations in which the α-bimodality arises. Using simulations from the NIHAO-UHD project, Buck (2020) found that the high-α sequence formed first, and the low-α sequence was formed after the interstellar medium was diluted with metal-poor gas brought in by a gas-rich merger. A chemical bimodality was also found in the VINTERGATAN simulation and linked to a rapid formation of an outer, metal-poor, low-|$\rm [\alpha /Fe]$| disc surrounding the inner, metal-rich, old high-|$\rm [\alpha /Fe]$| population (Agertz et al. 2021; Renaud et al. 2021a, b). Assuming two main starburst episodes ∼6 Gyr apart, Lian et al. (2020) found that the α-dichotomy can result from the rapid suppression of star formation after the first starburst. Conversely, Clarke et al. (2019) showed that an α-bimodality can occur even if the α-sequences form non-sequentially, where the high-α sequence forms from rapidly self-enriching star-forming clumps, while the low-α sequence is produced by radially distributed star formation. Khoperskov et al. (2021) also suggested that mergers are not necessary to create the bimodality, while Sharma et al. (2021) modelled the valley between the two sequences with time delays between SNII and SNIa. The latter interpretation, however, is at odds with the hybrid chemo-dynamical model by Minchev et al. (2013), where only a single gas infall was used and strong migration was present; yet, no separation into low- and high-α was found. Our results indicate that the high- and low-α sequences formed consecutively, and are consistent with work proposing that the low-α sequence was created by gas brought in by a gas-rich merger.

The fluctuations in the birth |$\rm \nabla [Fe/H](\tau)$| found at ∼4 and ∼6 Gyr (Fig. 2) create only small wiggles in the mono-|$\rm \text{R}_\text{birth}$| tracks in the [Fe/H]–[Mg/Fe] plane, and are able to reverse the continuous |$\rm [Fe/H]$| enrichment at each |$\rm \text{R}_\text{birth}$|⁠.

5.5 The weak possibility of chemical tagging in current-day data

Since the dynamics of stars evolve over time (Sellwood & Binney 2002; Roškar et al. 2008; Schönrich & Binney 2009b; Minchev & Famaey 2010), the goal of chemical tagging is to utilize the chemical homogeneity of stellar birth groups to determine stellar birth populations from the stars’ chemical abundances (Freeman & Bland-Hawthorn 2002; Bland-Hawthorn et al. 2010). With estimates of |$\rm \text{R}_\text{birth}$|⁠, we directly test this ability under current day observational data. We find that most |$\rm [X/Fe]$| distributions are normally distributed with different means for mono-age, mono-|$\rm \text{R}_\text{birth}$| populations, indicating that a given (age, |$\rm \text{R}_\text{birth}$|⁠) is reasonably chemically homogeneous (Fig. E1 in the Appendix). However, particularly for |$\rm \text{R}_\text{birth}$| close to the Galactic centre, there is a fair amount of overlap between the different age distributions. This means that under observational uncertainties, there are a variety of (age, |$\rm \text{R}_\text{birth}$|⁠) pairs for a given |$\rm [X/Fe]$|⁠. Using bins of 0.05 dex in |$\rm [Mg/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠, |$\rm [Mn/Fe]$|⁠, |$\rm [C/Fe]$|⁠, |$\rm [Al/Fe]$|⁠, and 0.1 dex in |$\rm [Ce/Fe]$| and |$\rm [Fe/H]$|⁠, we find there is an average age dispersion of 1 Gyr and an average |$\rm \text{R}_\text{birth}$| dispersion of 0.74 kpc for a given abundance pattern. Both of these values are on par with their measurement uncertainties, and suggest stronger versions of chemical tagging are infeasible until ages and abundances can be estimated more precisely (see also Casamiquela et al. 2021).

6 CONCLUSIONS

This work explored in detail the enrichment with lookback time of six chemical abundances across the Milky Way disc, using birth radii and stellar age for 87 426 APOGEE DR17 red giant stars. Since stars radially migrate away from their birth sites, a powerful tool in reconstructing the evolution of the Milky Way disc is to look at the preserved birth properties of stars – their abundances. With access to precise ages and |$\rm [Fe/H]$|⁠, |$\rm [X/Fe]$| for large samples of stars, we estimated the stellar birth radii without resorting to chemo-dynamical models, using the recent technique of Lu et al. (2022a). This allowed us to go beyond interpreting a diluted history from the present-day picture of the Galaxy and instead uncover the detailed chemical evolution of the Milky Way disc directly from the data.

In this work, we focused on representative elements from the α- (Mg, Ca), iron-peak (Mn), light (C), light odd-Z (Al), and s-process (Ce) families. While we performed quality cuts on signal-to-noise, abundance flags, and measurement uncertainties (Section 2), we remind the reader that the uncertainty in |$\rm [Ce/Fe]$| is larger than the uncertainty in the other abundances (0.1 dex compared to 0.03 dex) and is dependent on effective temperature (Cunha et al. 2017). Section  B in the Appendix briefly looks at the evolution of Ce using |$\rm [Ce/Fe]$| from the BAWLAS catalogue (Hayes et al. 2022), and the evolution of neutron capture elements are discussed in more detail in our work in preparation. Our key findings are summarized below:

• We found that the birth metallicity gradient evolution had three steepening phases ∼4, ∼6, and 9 Gyr ago (Fig. 2). This is in contrast to previous work using a similar method but finding monotonically flattening |$\rm \nabla [Fe/H](\tau)$| with time (Minchev et al. 2018) or only one fluctuation at ∼8 Gyr (Lu et al. 2022a). We tentatively associate these fluctuations in |$\rm \nabla [Fe/H](\tau)$| with dilution in |$\rm [Fe/H]$| from gas brought in by the GSE and Sagittarius mergers.

• The |$\rm \nabla [Fe/H](\tau)$| fluctuations at lookback times ∼4 and ∼6 Gyr result in flattening (and for some |$\rm \text{R}_\text{birth}$| and times even inversion) of the AMR for mono-|$\rm \text{R}_\text{birth}$| populations (Fig. 2, right). This breaks the expectation that a given |$\rm \text{R}_\text{birth}$| is monotonically increasing in |$\rm [Fe/H]$|⁠, even if no scatter is present at the time of birth (e.g. Minchev et al. 2018). See also fig. 8 in Buck (2020).

• There is structure in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| plane for the α-, light odd-Z, and iron-peak elements that is lost in the observable |$\rm \text{R}_\text{guide}$|–|$\rm [X/Fe]$| plane. This structure corresponds to the evolutionary events that create the high- and low-α sequences (Fig. 6).

• As a whole, the radial |$\rm [X/H]$| gradient flattens over time (Fig. 5). The shape of the |$\rm \text{R}_\text{birth}$|–|$\rm [X/H]$| radial gradient is most similar to the |$\rm \text{R}_\text{guide}$|–|$\rm [X/H]$| gradient for the youngest stars, while the older stars have significantly different gradients (Fig. 4). This is a consequence of radial mixing, as these younger stars have not had time to be strongly affected by radial migration.

• Except for |$\rm [Ce/Fe]$|⁠, the other |$\rm [X/Fe]$| and |$\rm [X/H]$| abundances examined in this work are good tracers of age once conditioned on |$\rm \text{R}_\text{birth}$|⁠. Mono-|$\rm \text{R}_\text{birth}$| populations do not monotonically increase in |$\rm [X/H]$| over lookback time due to the fluctuations in the metallicity gradient <8 Gyr ago (Fig. 9).

• We find that the gap in the |$\rm [Fe/H]$|–|$\rm [\alpha /Fe]$| plane comes from the quick transition between Type-II and Type-Ia Supernova dominated evolutionary phases that causes a fast drop in |$\rm [\alpha /Fe]$|⁠, in agreement with previous work (Lu et al. 2022a).

• In order to explain the decline in |$\rm [Mg/Fe]$| with increasing |$\rm [Fe/H]$| in the high-α sequence, one needs strong inside-out disc formation very early on (Fig. 7). This causes the transition from Type II to Type Ia Supernova in the inner few kpc sooner than in the outer radii, and can be linked to the strong positive [X/Fe] gradients of α-elements for the high-α sequence (top panels of Fig. 6).

Our results suggest that the Milky disc should be studied as one evolving entity, rather than splitting it into low- and high-α populations. We find that the α-sequences are not discrete discs, but rather they represent different times and birth radii, where the high-α sequence is simply composed of stars born at |$\rm \text{R}_\text{birth}$|<6 kpc with age >8 Gyr and the low-α sequence is a compilation of different ages and |$\rm \text{R}_\text{birth}$| (see Fig. 7). This conclusion is supported by the fact that except in chemical space, a bimodality in disc stars is not seen in kinematics or morphology (e.g. Bovy, Rix & Hogg 2012a).

The data-driven technique used in this work recovers the trends in chemical relations expected from detailed forward modelling (Matteucci & Francois 1989; Chiappini et al. 1997; Minchev et al. 2013; Buck 2020; Buck et al. 2020; Hemler et al. 2021; Johnson et al. 2021). Here, however, we infer these trends directly from the data using a completely different approach with less prior assumptions.

The detailed evolution of chemical abundances with birth radius and time that we derive in this work can serve as a prior at every time step for chemical evolution models, which are typically constrained by only the present-day metallicity gradient. Such work, however, needs to decipher the contribution from now-dead massive stars, galactic fountains, gas infall from cosmic filaments and satellites, etc., which, although all are naturally accounted for in our results, cannot be quantified by our approach.

ACKNOWLEDGEMENTS

Funding for the Sloan Digital Sky Survey V has been provided by the Alfred P. Sloan Foundation, the Heising-Simons Foundation, the National Science Foundation, 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 Carnegie Institution for Science, Chilean National Time Allocation Committee (CNTAC) ratified researchers, the Gotham Participation Group, Harvard University, Heidelberg University, The Johns Hopkins University, L’Ecole polytechnique fédérale de Lausanne (EPFL), Leibniz-Institut für Astrophysik Potsdam (AIP), Max-Planck-Institut für Astronomie (MPIA Heidelberg), Max-Planck-Institut für Extraterrestrische Physik (MPE), Nanjing University, National Astronomical Observatories of China (NAOC), New Mexico State University, The Ohio State University, Pennsylvania State University, Smithsonian Astrophysical Observatory, Space Telescope Science Institute (STScI), the Stellar Astrophysics Participation Group, Universidad Nacional Autónoma de México, University of Arizona, University of Colorado Boulder, University of Illinois at Urbana-Champaign, University of Toronto, University of Utah, University of Virginia, Yale University, and Yunnan University.

We acknowledge Lucy(Yuxi) Lu for helpful discussions. BR and IM acknowledge support by the Deutsche Forschungsgemeinschaft under the grant MI 2009/2-1. This work was partially funded by the Spanish MICIN/AEI/10.13039/501100011033 and by the ‘ERDF A way of making Europe’ funds by the European Union through grant RTI2018-095076-B-C21 and PID2021-122842OB-C21, and the Institute of Cosmos Sciences University of Barcelona (ICCUB, Unidad de Excelencia ‘María de Maeztu’) through grant CEX2019-000918-M. FA acknowledges financial support from MCIN/AEI/10.13039/501100011033 through grants IJC2019-04862-I and RYC2021-031638-I (the latter co-funded by European Union NextGenerationEU/PRTR).

DATA AVAILABILITY

The data sets used and analysed for this study are derived from data released by APOGEE DR17, Queiroz et al. (2023), and Anders et al. (2023). The rest of the relevant data sets are available from the corresponding author upon reasonable request.

Footnotes

1

In this study, we assume that most of the photospheric abundances reflect the abundance pattern in a star’s birth cloud.

2

We also examine Ce trends using |$\rm [Ce/Fe]$| provided in the BAWLAS catalogue. See Section  B in the Appendix.

3

Additionally, we note that the uncertainties of both the age (Anders et al. 2023) and the |$\rm \text{R}_\text{birth}$| determination methods (Section 3) grow with age, meaning that for older ages our method becomes less reliable.

4

It should be kept in mind that different age estimates are not necessarily on the same scale (depending on the method and the use of prior).

APPENDIX A: CONSISTENCY ACROSS AGE CATALOGUES

In this section, we show that the radial abundance gradients throughout time found in the main body of the text are not dependent on age derivation. The age catalogue of Anders et al. (2023) estimates age from stellar parameters and elemental abundance ratios using the supervised machine learning technique XGBoost. Here, we now use ages derived from isochrone fitting for APOGEE DR17 subgiant branch stars using the starhorse catalogue (Queiroz et al. 2023). We follow the same methods as in Section 3 to estimate |$\rm \text{R}_\text{birth}$| for 13 735 stars with |$|{\rm [Fe/H]}| \le 1$|⁠, |z| ≤ 1 kpc, eccentricity <0.5, |${\rm [Fe/H]}_\text{err} \lt 0.015$|⁠, and age_err < 1.25 Gyr.

After removing stars with flagged abundances and |${\rm [X/Fe]}_\text{err}\gt 0.05$| dex, we perform this comparison analysis for three elements – |$\rm [Mg/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠, and |$\rm [Mn/Fe]$| – of nearly 13 000 subgiant branch stars. Fig. A1 shows the |$\rm [X/Fe]$| radial trends as a function of both |$\rm \text{R}_\text{birth}$| and |$\rm \text{R}_\text{guide}$| for mono-age populations. Similar to Fig. 6, the inherent radial gradients tend to be steeper than looking at the present-day Galactic disc. The bimodal structure found in the |$\rm \text{R}_\text{birth}$|–|$\rm [\alpha /Fe]$| plane is also still preserved (Fig. A2). While the distribution of the low-α sequence in the |$\rm \text{R}_\text{birth}$|–|$\rm [Mg/Fe]$| plane contain minimal differences – possibly due to a difference in sample size or evolutionary state – the overall shape of this plane is quite similar. The similarities we find in all results – specifically the multiple fluctuations in |$\rm \nabla [Fe/H](\tau)$|⁠, the radial |$\rm [X/Fe]$| and |$\rm [X/H]$| gradients being steeper with |$\rm \text{R}_\text{birth}$| than |$\rm \text{R}_\text{guide}$|⁠, the structure in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| plane, mono-|$\rm \text{R}_\text{birth}$| populations lying on steeper tracks in the age–|$\rm [X/H]$| and age–|$\rm [X/Fe]$| planes, and the non-monotonic increase in [X/H] at a given |$\rm \text{R}_\text{birth}$| – show that the conclusions drawn in this work are not an artefact of stellar age catalogue or evolutionary state, and are indicative of the Milky Way’s chemical evolution.

Figure A1.

Comparison with isochrone ages derived with starhorse. Left 3 panels: Mono-age populations of 12 949 APOGEE DR17 subgiant branch stars in the Top: |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| and Bottom: |$\rm \text{R}_\text{guide}$|–|$\rm [X/Fe]$| planes. Right 3 panels: Mono-|$\rm \text{R}_\text{birth}$| populations in the Top: age–|$\rm [X/Fe]$| plane and Bottom: age–|$\rm [X/H]$| plane of the StarHorse sample. Consistencies seen between these planes and the ones in Figs 6 and 9 suggest that the conclusions drawn in Section 4 are not an artefact of stellar age derivation, and are representative of Milky Way evolution.

Figure A2.

The |$\rm \text{R}_\text{birth}$|–|$\rm [Mg/Fe]$| plane of Left: 87,426 APOGEE DR17 red giant stars with ages derived from abundances, and Right: 12,949 APOGEE DR17 subgiant branch stars with ages derived from isochrone fitting. Since StarHorse has few stars with ages <2 Gyr, the left figure has been cut to exclude the younger stars for a consistent comparison. Despite different evolutionary states and age derivations, the structure of the high- and low-α sequences is apparent in both figures, suggesting that the structure is real.

APPENDIX B: CERIUM RESULTS USING BAWLAS ABUNDANCES

In Section 4, we used |$\rm [Ce/Fe]$| abundances given by ASPCAP in the APOGEE DR17 catalogue. Here, we compare our results using |$\rm [Ce/Fe]$| provided in the BACCHUS Analysis of Weak Lines in APOGEE Spectra (BAWLAS; Hayes et al. 2022) catalogue. For this analysis, we keep disc (|z| ≤ 1 kpc, eccentricity <0.5) stars with |$-1\le \rm{\rm [Ce/Fe]}\le 0.5$|⁠, and a Ce empirical error cut <0.07 dex. The cut of |$\rm [Ce/Fe]$|<0.5 is chosen to remove s-process enhanced stars (Hayes et al. 2022) (though our results are robust to the inclusion of them). This leaves us with 36 239 stars for this analysis.

Fig. B1 provides the main figures reproduced with the new abundances. Overall, the results are consistent with the results shown in Section 4. In particular, |$\rm [Ce/Fe]$| has a weak radial birth gradient for older stars which becomes stronger with decreasing age (top left panel of Fig. B1). However, the existence of a break point is no longer noticeable for the youngest and oldest populations, but seems to still exist for the 3–10 Gyr populations. The [Ce/H] gradient also sees significant flattening with decreasing age (top middle panel of Fig. B1), in agreement with Fig. 4. In Fig. 9, we saw that the oldest end of each mono-|$\rm \text{R}_\text{birth}$| population in the age–|$\rm [Ce/Fe]$| plane decreased in |$\rm [Ce/Fe]$|⁠. Now, with |$\rm [Ce/Fe]$| derived using BAWLAS, we find even clearer separation of mono-|$\rm \text{R}_\text{birth}$| populations in the age–|$\rm [Ce/Fe]$| plane, which do not show this downward trend. However, the correlation between |$\rm [Ce/Fe]$| and age is still dependent on |$\rm \text{R}_\text{birth}$|⁠. The similarities between |$\rm [Ce/Fe]$| abundances in the APOGEE DR17 and BAWLAS catalogues suggest that the trends we find are not an artefact of pipeline, and are an inherent property of Ce evolution.

Figure B1.

Comparison of Ce results using |$\rm [Ce/Fe]$| abundances given in the BAWLAS catalogue. The running means of mono-age populations are in the top left: |$\rm \text{R}_\text{birth}$|–|$\rm [Ce/Fe]$|⁠, bottom left: |$\rm \text{R}_\text{guide}$|–|$\rm [Ce/Fe]$|⁠, top middle: |$\rm \text{R}_\text{birth}$|–[Ce/H], and bottom middle: |$\rm \text{R}_\text{guide}$|–[Ce/H] planes for 36 239 stars. The running means of mono-|$\rm \text{R}_\text{birth}$| populations in the age–|$\rm [Ce/Fe]$| and age–[Ce/H] planes are shown in the top right and bottom right panels, respectively. The overall trends are similar to Figs 4, 6, and 9, with mono-|$\rm \text{R}_\text{birth}$| populations showing more separation in the age–|$\rm [Ce/Fe]$| plane.

APPENDIX C: CONSISTENCY WHEN VARYING THE STRENGTH OF THE METALLICITY GRADIENT

As discussed in Section 3, we chose min(⁠|$\rm \nabla [Fe/H](\tau)$|⁠)  = −0.15 dex kpc⁻¹ based on previous works (Minchev et al. 2018; Lu et al. 2022a). The left panels of Fig. C1 show the distribution of |$\rm \text{R}_\text{birth}$| for mono-age populations (similar to Fig. 3) when |$\rm \text{R}_\text{birth}$| are derived with min(⁠|$\rm \nabla [Fe/H](\tau)$|⁠)  = −0.1, −0.15, and −0.2 dex kpc⁻¹. When the minimum |$\rm [Fe/H]$| gradient is −0.1 dex kpc⁻¹, the |$\rm \text{R}_\text{birth}$| distributions in the solar neighbourhood become too wide, and implies most of the 0–2 Gyr stars have migrated inwards, while a minimum |$\rm [Fe/H]$| gradient of −0.2 dex kpc⁻¹ indicates the younger aged stars have migrated outwards by over 3 kpc on average. Thus, our choice of min(⁠|$\rm \nabla [Fe/H](\tau)$|⁠)  = −0.15 dex kpc⁻¹ aligns most with expectations.

Figure C1.

The Left: distribution of |$\rm \text{R}_\text{birth}$| for stars in the solar neighbourhood, Middle: |$\rm [X/Fe]$| radial gradient over lookback time (Δlookback time  = 1 Gyr), and Right: |$\rm \text{R}_\text{birth}$|–|$\rm [Mg/Fe]$| plane when |$\rm \text{R}_\text{birth}$| is derived using a minimum metallicity gradient of Top row: −0.1 dex kpc⁻¹, Middle row: −0.15 dex kpc⁻¹ (used in the main body of this work), and Bottom row: −0.2 dex kpc⁻¹. The results presented in Section 4 are consistent for varying minimum |$\rm [Fe/H]$| radial gradients, with the overall trends in each plane being similar and the numerical values shifting.

The middle and right panels of Fig. C1 illustrate that even if our assumption of min(⁠|$\rm \nabla [Fe/H](\tau)$|⁠) is not optimal, the overall trends of our results are consistent when this assumption is varied. A more shallow gradient (min(⁠|$\rm \nabla [Fe/H](\tau)$|⁠)  = −0.1 dex kpc⁻¹) leads to flatter |$\rm [X/Fe]$| radial gradients over time, whereas a steeper gradient (min(⁠|$\rm \nabla [Fe/H](\tau)$|⁠)  = −0.2 dex kpc⁻¹) creates steeper |$\rm [X/Fe]$| radial gradients for a given lookback time. Despite the numerical differences between the |$\rm [X/Fe]$| radial gradients derived using different gradient strengths, the relative trends among the gradients is preserved. Similarly, the bimodal structure of the high- and low-α sequences in the |$\rm \text{R}_\text{birth}$|–|$\rm [Mg/Fe]$| plane found in Section 4.2 is still visible under different metallicity gradient strengths. The third clump of younger aged stars is more prominent under the steeper |$\rm \nabla [Fe/H](\tau)$|⁠, but can be seen using the different strengths.

APPENDIX D: RADIAL ABUNDANCE GRADIENTS WITH LOOKBACK TIME FOR MANY ABUNDANCES

Table D1 presents the numerical values of the |$\rm [X/Fe]$| radial gradients with lookback time. For each abundance, the number of stars used to determine the gradients is given in the second column. Only stars with unflagged |$\rm [X/Fe]$| with |$|{\rm [X/Fe]}| \lt 1$| and |${\rm [X/Fe]}_\text{err}\lt 0.1$| dex for |$\rm [Ce/Fe]$|⁠, or <0.05 dex for all other |$\rm [X/Fe]$| are used.

Table D1.

Numerical |$\rm [X/Fe]$| radial gradients over lookback time. Gradients are found by fitting a linear model to a mono-age population.

		|$\rm [X/Fe]$| radial birth gradient [dex kpc−1] for given lookback time [Gyr]
|$\rm [X/Fe]$|	n	0–1	1–2	2–3	3–4	4–5	5–6	6–7	7–8	8–9	9–10	10–11	11–12
[C/Fe]	144070	0.0	0.0	0.0	0.01	−0.0	−0.01	−0.01	−0.01	−0.0	0.01	0.01	−0.0
[N/Fe]	143986	−0.02	−0.02	−0.01	−0.02	−0.02	−0.03	−0.03	−0.02	−0.01	−0.01	−0.01	−0.01
[Al/Fe]	144072	0.0	0.01	0.01	0.02	0.02	0.02	0.02	0.03	0.04	0.04	0.03	0.02
[Mg/Fe]	144087	0.01	0.01	0.01	0.02	0.02	0.02	0.02	0.02	0.03	0.03	0.03	0.02
[O/Fe]	142406	0.01	0.01	0.02	0.03	0.02	0.02	0.03	0.03	0.04	0.05	0.04	0.04
[Si/Fe]	144090	0.0	0.0	0.01	0.01	0.01	0.01	0.01	0.02	0.02	0.03	0.03	0.02
[Ca/Fe]	144089	0.01	0.01	0.01	0.01	0.01	0.02	0.02	0.02	0.02	0.02	0.02	0.02
[Mn/Fe]	144090	−0.02	−0.02	−0.02	−0.02	−0.03	−0.03	−0.03	−0.03	−0.03	−0.03	−0.02	−0.02
[Cr/Fe]	89091	−0.0	−0.0	−0.01	−0.01	−0.01	−0.01	−0.01	−0.01	−0.01	−0.02	−0.02	−0.01
[Co/Fe]	33877	−0.04	−0.04	−0.03	−0.04	−0.04	−0.04	−0.03	−0.03	−0.02	−0.01	0.01	0.01
[Ce/Fe]	93123	0.06	0.04	0.04	0.04	0.03	0.04	0.03	0.01	0.01	−0.0	−0.02	−0.02

APPENDIX E: THE RELATIONSHIP BETWEEN AGE, BIRTH RADII, AND ABUNDANCES

Fig. 6 suggests that there is a relationship between |$\rm [X/Fe]$|⁠, age, and the |$\rm \text{R}_\text{birth}$| we derived in this work. In this section, we examine this relationship to test if abundances, age, and |$\rm \text{R}_\text{birth}$| have a unique link.

The left five columns of Fig.E1 show the |$\rm [X/Fe]$| distribution conditioned on age and birth radius alongside the |$\rm [X/Fe]$| distribution for the entire sample (black dashed line). Even when the stellar sample shows skew (⁠|$\rm [Al/Fe]$|⁠, |$\rm [Mn/Fe]$|⁠), long tails (⁠|$\rm [Mg/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠), or other non-normal traits (⁠|$\rm [C/Fe]$|⁠), most |$\rm [X/Fe]$| for a given age and |$\rm \text{R}_\text{birth}$| are normally distributed with little variance. This implies that a given birth time and place is reasonably chemically homogeneous in these abundances. The few non-normally distributed |$\rm [X/Fe]$| (e.g. long tails in |$\rm [Mg/Fe]$| for |$2 \lt \rm{\rm \text{R}_\text{birth}}\lt 4$| kpc, bimodality in |$\rm [Al/Fe]$| for |$4 \lt \rm{\rm \text{R}_\text{birth}}\lt 6$| kpc) tend to be for older aged stars. These traits may be a consequence of poor age estimates of older stars (as it is difficult to measure the difference between an e.g. 10 and 12 Gyr old star; Mackereth et al. 2019), or that period in time corresponds to a quickly enriching environment in the inner disc. We note that the 8–10 Gyr age bin encompasses both high- and low-α stars.

Figure E1.

Left five columns: |$\rm [X/Fe]$| distributions for mono-age, mono-|$\rm \text{R}_\text{birth}$| populations alongside the overall |$\rm [X/Fe]$| distribution for the entire sample (black dashed line). Right: |$\rm [X/Fe]$| distributions for mono-age populations alongside the distribution for the entire sample (black dashed line). Most abundances show minimal separation when only conditioning on age (e.g.|$\rm [Mn/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠, |$\rm [Al/Fe]$|⁠), but show distinct differences when additionally conditioning on place of birth.

The old (10–12 Gyr) stars show a clear separation from the other age distributions for all mono-|$\rm \text{R}_\text{birth}$| populations for |$\rm [Al/Fe]$|⁠, |$\rm [Mg/Fe]$|⁠, |$\rm [Mn/Fe]$|⁠, and |$\rm [Ca/Fe]$|⁠, and across larger |$\rm \text{R}_\text{birth}$| for |$\rm [C/Fe]$|⁠. This distinction roughly corresponds to the high- and low-α sequences, and demonstrates that the sequences are most similar in |$\rm [X/Fe]$| at smaller |$\rm \text{R}_\text{birth}$|⁠. This idea is clear in the |$\rm [Fe/H]$|–|$\rm [Mg/Fe]$| plane, where the sequences converge for small |$\rm \text{R}_\text{birth}$| and have no obvious distinction (Fig. 7). For larger |$\rm \text{R}_\text{birth}$|⁠, some elements have an average separation between the older and younger age groups of 0.2 dex, with minimal stars between the oldest and younger age groups (e.g. |$\rm [Mg/Fe]$| and |$\rm [Ca/Fe]$| at |$6\lt \rm{\rm \text{R}_\text{birth}}\lt 8$|⁠). Fig. E1 clearly illustrates that the older age stellar population has a different chemical composition than the younger populations, which suggests that there was a major enrichment event that quickly diluted the interstellar medium and caused a difference in these abundances.

To demonstrate the additional information gained with knowledge of |$\rm \text{R}_\text{birth}$|⁠, the rightmost column of Fig. E1 presents each |$\rm [X/Fe]$| distribution only conditioned on stellar age. For half of the abundances considered in this work (⁠|$\rm [Al/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠, and |$\rm [Mn/Fe]$|⁠), there is very little separation between age bins beyond a bulk difference between the oldest (10–12 Gyr) and younger aged populations. For these elements, the <10 Gyr age bins peak at the same location as the overall population (black dashed line), with only the 10–12 Gyr age bin having a wider, sometimes bimodal, distribution with a different peak. |$\rm [Mg/Fe]$| also has a similar trend to these three abundances, though the younger aged populations have slight variations in the location of their peak. Once utilizing birth radii information, the mono-age populations separate as a function of |$\rm \text{R}_\text{birth}$|⁠. The |$\rm [Mg/Fe]$| bimodality of the 10–12 Gyr population becomes, for most |$\rm \text{R}_\text{birth}$|⁠, more normally distributed with less variation. The younger aged populations also begin to show less variability in |$\rm [Mg/Fe]$| when additionally conditioning on |$\rm \text{R}_\text{birth}$|⁠, with the peak of each mono-age, mono-|$\rm \text{R}_\text{birth}$| population being slightly shifted.

Opposite to the four abundances listed so far, the light element abundance |$\rm [C/Fe]$| behaves differently than the iron-peak and α-elements. With only conditioning on age, |$\rm [C/Fe]$| shows clear separation between different age bins, and additionally conditioning on |$\rm \text{R}_\text{birth}$| does not reveal any additional distinction between the groups. This again highlights the strong relationship between |$\rm [C/Fe]$| and age. The s-process abundance |$\rm [Ce/Fe]$| also behaves differently than the iron-peak and α-elements. |$\rm [Ce/Fe]$| has very little difference in the older mono-age populations, with the younger (<4 Gyr) stars showing the most deviation. Once binned in |$\rm \text{R}_\text{birth}$|⁠, the older populations are separated for smaller |$\rm \text{R}_\text{birth}$|⁠, and show minimal separation for |$\rm{\rm \text{R}_\text{birth}}\gt 2$| kpc. All the mono-age, mono-|$\rm \text{R}_\text{birth}$| populations show no differences in |$\rm [Ce/Fe]$| for |$4 \lt \rm{\rm \text{R}_\text{birth}}\lt 6$| kpc. This location in the disc is where we find the pivot point of the |$\rm \text{R}_\text{birth}$|–|$\rm [Ce/Fe]$| radial gradients (Fig. 6), suggesting that the inner and outer disc had different |$\rm [Ce/Fe]$| enrichment histories. The second row of Fig. E1 also depicts that the lower |$\rm [Ce/Fe]$| values are in older stars produced in the inner disc, while the higher |$\rm [Ce/Fe]$| values are in younger stars born in the outer disc.

Within each mono-|$\rm \text{R}_\text{birth}$| population, the relative distribution of the mono-age populations reveals the |$\rm [X/Fe]$| time gradient at a given |$\rm \text{R}_\text{birth}$|⁠. As discussed above, the mono-age populations show distinction in the light element |$\rm [C/Fe]$|⁠, with minimal additional separation when also conditioning on place of birth. Within a mono-|$\rm \text{R}_\text{birth}$| population, the |$\rm [X/Fe]$| time gradient is fairly consistent, particularly for |$\rm{\rm \text{R}_\text{birth}}\gt 6$| kpc. For |$\rm{\rm \text{R}_\text{birth}}\lt 6$| kpc, the younger populations have higher |$\rm [C/Fe]$|⁠, showing a weaker gradient with time.

|$\rm [Al/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠, and |$\rm [Mg/Fe]$| also show weaker time gradients closer to the Galactic centre, with the gradients steepening at larger |$\rm \text{R}_\text{birth}$|⁠. Specifically, for |$\rm{\rm \text{R}_\text{birth}}\lt 2$| kpc, the mono-age populations are clustered about the peak of the non-conditional sample (black dashed line), showing that the |$\rm [X/Fe]$| distribution does not vary much near the Galactic centre. For larger |$\rm \text{R}_\text{birth}$|⁠, the difference between |$\rm [X/Fe]$| for different aged stars becomes clearer, especially between the 10–12 Gyr age group and the others. Even within the <8 Gyr populations, the distributions peak at slightly different locations, indicating a smoothly continuous evolution with time. We also find for these three abundances that for any given |$\rm \text{R}_\text{birth}$|⁠, a wide range of |$\rm [X/Fe]$| values can be created.

While the α- and light elements show weaker time gradients for smaller |$\rm \text{R}_\text{birth}$|⁠, |$\rm [Mn/Fe]$| does the opposite. For |$\rm{\rm \text{R}_\text{birth}}\lt 4$| kpc, the different age populations are structured with time. However for |$\rm{\rm \text{R}_\text{birth}}\gt 6$| kpc, the <8 Gyr stellar groups show little variation with lookback time. This indicates that nearly all values of |$\rm [Mn/Fe]$| are formed at lower radii, but only stars with a small range of |$\rm [Mn/Fe]$| are born further out in the disc.

APPENDIX F: EXTRA PLOTS

In this section, we provide additional plots. Fig. F1 indicates the cut we use to define the high- and low-α sequences by eye. Fig. F2 illustrates the robustness of the fluctuations in |$\rm \nabla [Fe/H](\tau)$| by showing the gradient recovered using GALAH DR3 (Buder et al. 2021) subgiant branch disc stars with ages derived by isochrone fitting from the starhorse catalogue (Queiroz et al. 2023). Fig. F3 shows the linear fits for each mono-age population in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| plane, with |$\rm [C/Fe]$|⁠, |$\rm [Ce/Fe]$|⁠, and |$\rm [Mn/Fe]$| being fit with two linear models (one before and one after the defined break points).

Figure F1.

2D density distribution of 87 426 APOGEE DR17 red giant stars in the [Fe/H]–[Mg/Fe] abundance plane. The dashed teal line represents our simple by eye division into the high- and low-α sequence.

Figure F2.

Inferred metallicity gradient (red) and scatter in |$\rm [Fe/H]$| (blue) using GALAH DR3 subgiant branch disc (|z| < 1, eccentricity <0.5, |$|{\rm [Fe/H]}| \lt 1$|⁠) stars with ages derived by isochrone fitting from the starhorse catalogue. We additionally remove flagged |$\rm [Fe/H]$| and age measurements and keep stars with |$\rm [Fe/H]$|_err < 0.1 and age₈₄ − age₁₆ < 3. We calculate the gradient using a bin size of 1 Gyr every 0.2 Gyr (similar as done in Fig. 10). The three fluctuations found in Section 3 are again recovered, showing their robustness across survey, age derivation, and evolutionary state.

Figure F3.

Running mean of mono-age populations in the |$\rm \text{R}_\text{birth}$|–|$\rm [X/Fe]$| plane (solid lines) with the estimated fits (dashed lines). Mono-aged populations are fit as a linear function in |$\rm [Mg/Fe]$|⁠, |$\rm [Ca/Fe]$|⁠, and |$\rm [Al/Fe]$|⁠, while the mono-age groups show non-linear trends between |$\rm \text{R}_\text{birth}$| and |$\rm [C/Fe]$|⁠, |$\rm [Ce/Fe]$|⁠, and |$\rm [Mn/Fe]$|⁠. Therefore we fit these three abundances as a function of two linear functions, one before a break point and one after. The break points are different for each mono-age group, as well as each |$\rm [X/Fe]$|⁠.