© ESO 2019

1 Introduction

The combination of precise chemical abundances with stellar ages in different stellar populations opens the door to a more insightful view on Galactic archaeology. Moreover, there is an active discussion on how the thin and thick disk components of the Galaxy should be defined, with an increasing number of works suggesting that age rather than kinematics is a better parameter to differentiate both populations (e.g., Haywood et al. 2013; Bensby et al. 2014). Therefore, a great effort to derive reliable ages in big samples of stars has been done in the previous years. The use of asteroseismic observations represents a significant advance in the derivation of accurate ages but the samples analyzed with this method are still limited in size (e.g., Kepler LEGACY stars in Nissen et al. 2017) or limited to a range of stellar parameters (e.g., red giants in APOKASC and in CoRoGEE dataset, Pinsonneault et al. 2018; Anders et al. 2017, respectively) in order to allow for a global comprehensive analysis of the Milky Way. However, the arrival of Gaia data will have a significant impact in the derivation of reliable ages for large sample of stars across the full Galaxy. Indeed, some large spectroscopic surveys are already taking advantage of Gaia DR1 and DR2 to derive ages and distances, such as LAMOST (e.g., Tian et al. 2018; Yu & Liu 2018), Gaia-ESO survey (Randich et al. 2018), GALAH (e.g., Buder et al. 2019) and APOGEE

Furthermore, constraining the temporal evolution of different chemical species can help to understand which nucleosynthesis channels aretaking place at different ages and evaluate the relative importance of stellar yields at a given time. Moreover, the different sources producing elements (e.g., massive stars, neutron star mergers, AGB, SNeIa) are not distributed in an homogeneous way across the Galaxy. Therefore, the current abundances of a given stellar population are the outcome of several conditions that need to be accounted for to reconstruct the chemical evolution of the Galaxy. As shown in this work, the distinction of different abundance-age trends requires detailed and very precise chemical abundances that are not always possible to achieve in larger spectroscopic surveys. Therefore, although our sample is of a modest size compared to previously mentioned surveys and only covers the solar neighborhood, the high resolution and high signal-to-noise ratio (S/N) of the data will serve to make an important contribution for future models of Galactic chemical evolution (GCE).

In this work, we provide stellar ages for the HARPS-Guaranteed Time Observations (GTO) sample in order to study the temporal evolution of chemical species with different nucleosynthetic origin and to analyze the feasibility of using different abundance ratios to estimate stellar ages. The structure of the paper is the following: in Sects. 2 and 3 we describe the sample of chemical abundances and the derivation of stellar ages. In Sects. 4 and 5 we evaluate the temporal evolution of metallicity and different abundances ratios. Section 6 is devoted to find statistically significant correlations between age and abundance ratios with the aim of obtaining simple relations to derive stellar ages. Finally, in Sect. 7 we present the summary of the results.

2 Stellar parameters and chemical abundances

The baseline sample used in this work consist of 1111 FGK stars observed within the context of the HARPS-GTO planet search programs (Mayor et al. 2003; Lo Curto et al. 2010; Santos et al. 2011). The final spectra have a resolution of R ~115 000 and high S/N (45% of the spectra have 100 <S/N < 300, 40% of the spectra have S/N > 300 and the mean S/N is 380).

Precise stellar parameters for the full sample of 1111 stars within the HARPS-GTO program were homogeneously derived in Sousa et al. (2008, 2011a,b). The parameters for cool stars were revised by Tsantaki et al. (2013) using a special list of iron lines which was later applied to the full sample in Delgado Mena et al. (2017; hereafter DM17), also correcting the spectroscopic gravities. From the 1111 stars in the original sample, the derivation of parameters converged to a solution for 1059 of them. Our stars have typical T_eff values between 4500 and 6500 K and surface gravities mostly lie in the range 4 < logg < 5 dex meanwhile the metallicity covers the region −1.39 < [Fe/H] < 0.55 dex. Chemical abundances of Cu, Zn, Sr, Y, Zr, Ba, Ce, Nd, and Eu were determined under local thermodynamic equilibrium (LTE) using the 2014 version of the code MOOG (Sneden 1973) and a grid of Kurucz ATLAS9 atmospheres (Kurucz 1993). For more details about the sample and the analysis we refer the reader to DM17. In that work we also provide updated values of Mg, Al, Si, Ca, and Ti using the corrected values of log g and T_eff and the EWs measured by Adibekyan et al. (2012).

In the analysis presented in DM17 we used as solar reference abundances those derived with the spectroscopic log g (4.43 dex) instead of using the corrected value considering its T_eff (4.39 dex, obtained with Eq. (2) of such work) as done for the full sample. The differences are minimal and well below the errors (<0.006 dex for neutral elements and <0.018 dex for ionized species) so they do not change the conclusions presented in DM17. However, since our aim is to obtain empirical relations to derive stellar ages through multivariable linear regressions we must consider this small difference. The chemical clocks are typicallymade by the subtraction of one neutral element from one ionized element (e.g., [Y/Mg]) so the difference might be slightly bigger and could produce an small offset with the empirical relations found by other authors (P. E. Nissen, priv. comm.). Therefore, we re-derived the [X/Fe] of the full sample by using the solar reference values obtained with the corrected log g. In Table 1 the solar abundances from DM17 are compared with the new ones using log g = 4.39 dex.

3 Derivation of stellar ages

In this work, we have derived three sets of masses and ages, the last two for comparison purposes:

• Gaia DR2parallaxes together with PARSEC isochrones;

• HIPPARCOS parallaxes together with PARSEC isochrones;

• T_eff and log g together with Yonsei–Yale isochrones.

For thefirst two sets we derived the masses, radii and ages with the PARAM v1.3 tool¹ using the PARSEC isochrones (Bressan et al. 2012) and a Bayesian estimation method (da Silva et al. 2006) together with the values for T_eff and [Fe/H] from DM17, the V magnitudes from the main HIPPARCOS catalog (Perryman et al. 1997) and the parallaxes from the second release (DR2) of Gaia (Gaia Collaboration et al. 2016, 2018; Lindegren et al. 2018), available for 1057 out of 1059 stars. All the starsin our sample have errors in parallax well below the value of the parallax itself, with most of them having an error in parallax lower than 0.1 mas. The Bayesian inference is applied taking into account priors for the initial mass function (Chabrier 2001) and a constant Star Formation Rate. We used as a prior a maximum age of 13.5 Gyr. No correction for interstellar reddening was considered as all stars are in close distance. We note that the errors in T_eff and [Fe/H] that need to be input in PARAM to derive the ages should be absolute errors. Therefore, the final errors are the quadratic sum of the precision errors (those reported in DM17) and a systematic error of 60 K and 0.04 dex as determined by Sousa et al. (2011a).

The distribution of ages in the full sample can be observed in the upper panel of Fig. 1 for the different populations as defined in Adibekyan et al. (2011, 2013). We recall that the thin and thick disk stars are defined based in the chemical separation in [α/Fe] (being α the average of Si, Mg, and Ti) across different metallicities bins (see Figs. 1 and 2 and Fig. 9 in Adibekyan et al. 2011, and DM17, respectively). That work also revealed the existence of a high-α metal-rich (hereafter hαmr) population, with [Fe/H] >−0.2 dex and enhanced [α/Fe] ratios with respect to the thin disk. On the other hand, halo stars are defined based on their kinematics alone. Thick disk stars have a peak at around 12 Gyr and hαmr stars are clearly older on average than thin disk stars, with two peaks around 6 and 9 Gyr. As expected, most of the halo stars are older than 11 Gyr. We note that if we consider those stars with small errors in age (<1.5 Gyr) the distribution is different (see lower panel of Fig. 1), with thin disk stars peaking at younger ages because hotter stars tend to have smaller errors in parameters and thus also smaller errors in age. Also, the peak of thin disk and hαmr stars around 6 Gyr disappears since that was produced by cool stars with large errors in age. Finally, the majority of thick disk and halo stars with ages lower than ~8 Gyr are also removed due to their large errors in age.

As a comparison, we also derived masses and ages using parallaxes from the HIPPARCOS new reduction (van Leeuwen 2007) with the same aforementioned method. HIPPARCOS provides parallaxes for 1051 out of the 1059 stars within our sample. The third set of ages and masses was derived by using the q² Python package² (Ramírez et al. 2014) that considers Yonsei–Yale isochrones (Kim et al. 2002) and spectroscopic parameters (T_eff, log g and [Fe/H]). We note that for this set of ages we have used the corrected values of log g presented in DM17. In Fig. 2 we compare the results we obtain with the three previously mentioned methods. The masses obtained with different methods and/or parallaxes agree very well and the differences are quite small (⟨M_q2 – M_Gaia⟩ = 0.003 ± 0.051 M_⊙ and ⟨M_HIP – M_Gaia⟩ = –0.004 ± 0.030 M_⊙) although the comparison between Gaia masses and q² masses show some oscillations around the mean.

However, the situation is different for the ages for which we can find very large differences for some stars although the average differences are not very large: (⟨Age_q2 − Age_Gaia⟩ = 1.12 ± 1.96 Gyr and ⟨Age_HIP – Age_Gaia⟩ = –0.46 ± 1.28 Gyr). For example, there are many stars clustered around 8–10 Gyr in the q² ages dataset. They are very cool dwarfs and have similar stellar parameters (T_eff, log g and [Fe/H]) so the code delivers very similar ages. However, when using parallaxes and magnitudes the degeneracy is broken and their age range increases. Nevertheless, although the differences in age between the three methods present a moderate dispersion there are no large systematic differences. It is also clear that there is a large group of stars whose ages from Gaia are much larger than with HIPPARCOS or q² (see upper panel of Fig. 2). Most of the stars that have an age with Gaia greater than with HIPPARCOS also have smaller parallaxes (see Fig. 3) and the difference in age is not correlated with the magnitude.

Several works have determined the systematic offset between Gaia DR2 parallaxes and other samples with independently derived parallaxes, all of them finding that Gaia parallaxes are lower by 0.03–0.09 mas. However, such offset becomes more important for more distant stars. For example, in Stassun & Torres (2018) the comparison between eclipsing binaries parallaxes and Gaia DR2 shows that they are almost the same down to 4 mas. Since in our sample only three stars have a parallax smaller than 4 mas we decided to add a conservative value of 0.03 mas as suggested by the Gaia collaboration (Lindegren et al. 2018). Moreover, we also increased the errors in parallaxes to consider the ~30% underestimation in uncertainties for bright stars (Luri et al. 2018; Arenou et al. 2018). However, we remark that the increase of errors by such amount hardly affects the derived ages. We note that the average parallax difference for the 1047 stars in common between Gaia DR2 and HIPPARCOS is –0.239 ± 1.206 mas, greater than the comparison of the full sample made by Arenou et al. (2018; –0.118 ± 0.003 mas). In any case, we decided to use ages derived with Gaia DR2 parallaxes as final ages (given that those parallaxes are much more precise than HIPPARCOS). In order to have more reliable results we will consider in next subsections the 354 stars with errors in age lower than 1.5 Gyr (265 thin disk stars, 15 hαmr stars, 70 thick disk stars and 4 halo stars). These 354 stars have a large range in parameters T_eff : 5010−7212 K (95% between 5300 and 6500 K), logg: 3.73−4.71 dex (92% between 4.0 and 4.6 dex), [Fe/H]: –1.15−0.55 dex (93% between –0.7 and 0.4 dex). This chosen limit in age error is a compromise between having a reliable set of stellar ages but still with a sufficient number of stars for a meaningful analysis. However, we note that by cutting the age error at 1 or 2 Gyr instead of 1.5 Gyr the final conclusions of this work would not change. The final ages and masses are listed in Table 2.

Fig. 1 Normalized distribution (per stellar population) of ages for the full sample (upper panel) and a subset with small errors in age (lower panel) using Gaia DR2 parallaxes. The yellow, red, blue and purple histograms represent the thin, hαmr, thick and halo stars, respectively.

Fig. 2 Comparison of ages and masses obtained using Gaia DR2 or HIPPARCOS parallaxes and obtained with package q2. The red dots are the stars with errors in age lower than 1.5 Gyr. The green dashed lines shows the zero differences and the blue lines show the mean and standard deviation.

Fig. 3 Comparison of age differences vs. parallax differences between Gaia DR2 and HIPPARCOS. The V magnitudes are shown in a color scale.

4 Temporal evolution of [Fe/H] and [α/Fe]

In the upperpanel of Fig. 4 we depict the metallicity-age relation in the full sample. There is an overdensity of stars around 5–6 Gyr caused by most of the cool stars (≲5100 K) for which the errors in age are typically above 4 Gyr. Therefore, we show in the lower panel of the same figure the subsample with more precise ages (error in age lower than 1.5 Gyr). It is clear from the figure the large dispersion of ages at a given metallicity as previously observed in the solar neighborhood (e.g., Haywood et al. 2013; Bensby et al. 2014; Bergemann et al. 2014) or in other specific samples such as the dwarf stars in GALAH (Buder et al. 2019)or the giant stars in Kepler (Silva Aguirre et al. 2018), LAMOST-Kepler (Wu et al. 2018) or APOGEE (Feuillet et al. 2018). This confirms the weak age-[Fe/H] correlation first pointedout by Edvardsson et al. (1993) which is assumed to be caused by radial migration (Sellwood & Binney 2002). In the absence of radial migration we would expect to have a stronger age-[Fe/H] correlation. Furthermore, the dispersion in metallicity increases with age being 0.13 dex for stars younger than 2 Gyr, 0.19 dex between 2 and 4 Gyr and around 0.26 dex for older stars. However, there is a decrease of the dispersion to 0.16 dex for stars older than 12 Gyr (i.e., thick disk stars) which might be caused by the low number of metal-poor stars in our sample. Recent observations seem to suggest no strong radial metallicity gradient for the thick disk (Cheng et al. 2012; Hayden et al. 2015; Allende Prieto et al. 2006). This would mean that either there was not a radial metallicity gradient at the formation of the thick disk (meaning that the stars were formed from a well mixed material) or that there was a gradient that was flattened out by efficient radial migration. However, efficient radial migration in the thick disk would make the metallicity dispersion wider (e.g., Minchev et al. 2013) which, if real, contradicts our observed trend. Therefore, our observations together with a current flat metallicity gradient suggests that the thick disk was probably formed from a well mixed material (Haywood et al. 2015) and with no radial metallicity gradient.

The average [Fe/H] per age bin keeps around 0.08 dex up to ~3 Gyr, then it decreases to solar value for stars up to ~6 Gyr. At ages ≳6 Gyr, however, there is a more clear decreasing trend of [Fe/H] with age, reaching average values of −0.57 dex (12–13 Gyr) and −0.75 dex (>13 Gyr). The most metallic star (HD 108063, [Fe/H] = 0.55) has an age of 1.81 ± 0.06 Gyr but the majority of most metal-rich (with [Fe/H] ~ 0.4 dex) stars are indeed close to 4 Gyr. It is very probable that the most metal-rich stars were formed in the inner Galaxy and because they are older, they had time to radially migrate into the solar neighborhood (e.g., Sellwood & Binney 2002). Indeed, a recent study by Minchev et al. (2018) using our sample shows that the most metal rich stars have the smallest birth radii.

If we consider only those stars with a low error on age (<1.5 Gyr; middle panel of Fig. 4) we can see that only a few thick disk stars have ages lower than 11 Gyr and hαmr stars are always older than 4 Gyr. On the other hand very few thin disk stars are older than 10 Gyr in general agreement with previous works (e.g., Haywood et al. 2013) but in contrast with the recent works by Wu et al. (2018) and Silva Aguirre et al. (2018), where ages have been derived using asteroseismic constraints. We can also see a certain pile up of stars around 13.5 Gyr. This is an artifact of the prior we input as maximum age to avoid having ages older than the age of the Universe.

In the lower panel of Fig. 4 the relation of [O/H] is shown. We have taken the oxygen abundances from Bertran de Lis et al. (2015), derived for the HARPS-GTO sample. We note that those abundances were calculated with the previous set of parameters so we have rederived them here using the corrected log g from DM17³ although we note that the changes are not substantial. Oxygen is considered a “pure” α element, mostly produced by massive stars and ejected in SNeII (e.g., Woosley & Weaver 1995; Andrews et al. 2017) and its production along the time is quite constant in comparison with iron, with a maximum dispersion of 0.2 dex, around 6 Gyr, and average [O/H] values between –0.1 and 0.1 dex across all the ages.

Several works have proposed that [α/Fe] is a better indicator of age than [Fe/H] and this is also confirmed by our results (see the upper panel of Fig. 5 were α is the average of Mg, Si, and Ti). In contrast with metallicity, the dispersion of [α/Fe] ratios is much lower, with values of 0.03 dex up to 4 Gyr, 0.04 dex from 4 to 8 Gyr, 0.06 dex from 8 to 12 Gyr and 0.04 dex for stars older than 12 Gyr. The abundances of α elements with respect to iron in thin disk and hαmr stars show a clear increasing linear trend with age (with Spearman correlation coefficients, ρ, of 0.78 and 0.72, respectively). The tail of the distribution at old ages, formed by thick disk stars, also presents an increasing trend, though with a somewhat lower ρ of 0.59. This lower significance is mostly produced by the shorter age range of our chemically defined thick disk stars. In other works there is not distinction between thick disk and hαmr stars, thus having thick disk stars as young as 8–9 Gyr (e.g., Haywood et al. 2013; Bensby et al. 2014). If we consider the hαmr stars as part of the thick disk, the Spearman correlation coefficient for the [α/Fe]-age relation would be 0.74. Moreover, the lower significance of the correlation is also partly caused by the existence of a few “young” stars (~5–10 Gyr) with a similar α content as their older counterparts. In contrast, this region of intermediate ages and high-α abundances is well populated in the works by Wu et al. (2018) and Silva Aguirre et al. (2018) producing a rather flat behavior with age for the high-α sequence. Several works have found the existence of young α-rich stars (e.g., Fuhrmann et al. 2012; Fuhrmann & Chini 2017; Haywood et al. 2015), which seem to be more common in the inner disk of the Galaxy (Chiappini et al. 2015). These objects might be evolved blue stragglers (e.g., Fuhrmann et al. 2011; Jofré et al. 2016) or have been formed near the end of the Galactic bar (Chiappini et al. 2015).

Regarding the shape of the [α/Fe]-age trend, the work by Haywood et al. (2013)-that uses our HARPS-GTO sample of abundances-finds a change of slope around 8 Gyr, at which older stars show a steeper trend of [α/Fe] with age. In that work the oldest thin disk stars are about 10 Gyr, similar to the present study. However, in the Kepler and LAMOST-Kepler samples there are a significant number of older thin disk stars that follow the trend of younger stars and thus avoid the change of slope found in Haywood et al. (2013) and also found in the present study but at an older age. Nevertheless, we note here that the trends found by Wu et al. (2018) and Silva Aguirre et al. (2018) might be difficult to compare with those in the present study or in the work by Haywood et al. (2013) due to the much larger errors in age for old stars in those two samples of red giants.

In Fig. 6 we show the same relations for [α/Fe] using the ages from q², that are derived with the same isochrones as in Haywood et al. (2013) and making use of log g instead of parallaxes⁴. The trends are quite similar as those using ages derived with Gaia parallaxes. However, there are hardly no thin disk stars above 10 Gyr and [α/Fe] increases more steeply with age after ~8 Gyr whereas when using Gaia ages the change of slope seems to happen around 10–11 Gyr. In the dwarfs sample within GALAH presented by Buder et al. (2019) there are also many stars older than 11 Gyr with low [α/Fe] as in the works by Wu et al. (2018) and Silva Aguirre et al. (2018). Nevertheless, those ages might not be very precise due their large errors as reported by Buder et al. (2019). However, the [α/Fe] ratios by Buder et al. (2019) seem to present a steep trend with age for older ages (still with a high dispersion) in contrast with the work by Wu et al. (2018) and Silva Aguirre et al. (2018). Whether the steeper increase of [α/Fe] with age for older stars⁵ is taking place only in the solar neighborhood or if it is extendable to other parts of the Galaxy will be better understood with future releases of large spectroscopic surveys. In addition, the use of different methods to derive ages (e.g., asteroseismic constrains, using log g instead of parallax, applying a Bayesian approach instead of a χ² fitting) or the use of different sets of isochrones might produce different trends as shown in Figs. 5 and 6. The works by Haywood et al. (2013) and Bensby et al. (2014) both use Yonsei–Yale isochrones and they find the change of slope in the [α/Fe]-age relation around 8–9 Gyr in a similar way as our sample with q² ages (which also uses Yonsei–Yale). However, our ages obtained with PARSEC isochrones tend to be older which might be the reason for the change of slope occurring ~2 Gyr later and the fact of having a handful of thin disk stars older than 10 Gyr. In a similar way, the works by Anders et al. (2018), Minchev et al. (2018) that derive ages for our sample using the PARSEC isochrones within the StarHorse package, find a change of the [α/Fe]-age slope at ~10 Gyr. On the other hand, the work by Buder et al. (2019) uses the Dartmouth isochrones which also produce larger ages than Yonsei-Yale isochrones (see Fig. 3 in Haywood et al. 2013). This might explain why that work finds so many thin disk stars older than 11 Gyr and why the change of slope of the [α/Fe]-age seems to happen at a bit older age than that of Bensby et al. (2014; although still compatible within the uncertainties, see Fig. 14 in Buder et al. 2019).

The correlation between [α/Fe] and [Fe/H] is depicted in the lower panel of Fig. 5 with a color scale for the ages. Most of the stars with [α/Fe] ≳ 0.2 dex (92%) have ages greater than 10 Gyr and belong to the thick disk. This fact shows that the formation of the thick disk was fast and took place mostly before the thin disk. On the other hand, the thin disk stars show a great dispersion in ages, which increase as we move to lower [Fe/H] when considering stars of a similar [α/Fe]. If we look at the oldest thin disk stars, they have a similar age and [α/Fe] (around 0.1 dex) as the oldest hαmr (orange–green squares) but they are more metal-poor. Haywood et al. (2013) proposed that this group of older thin disk stars must have been formed in a different part (specifically the outer disk) than the hαmr stars⁶ since they have different metallicities and a higher rotational component in their velocities. The study of the AMBRE:HARPS sample by Minchev et al. (2018) also shows that those metal-poor thin disk stars have the largest birth radii. We can also see that there is a very clear stratification in ages (older as [α/Fe] increases) when looking at stars at a given [Fe/H], running from <2 Gyr for the most α-poor (at [Fe/H] ≳ –0.2 dex) to 8–10 Gyr as we reach the hαmr stars. Therefore, hαmr stars are clearly separated from the metal-rich thin disk both in terms of [α/Fe] and age. A similar conclusion is also reached with the analysis by Anders et al. (2018). Since most of our thick disk stars have ages above 12 Gyr it is difficult to see a temporal evolution in the thick disk as metallicity diminishes. However, in the lower panel of Fig. 6 where the same plot is done using the q² ages, thick disk stars tend to be older as [α/Fe] increases and [Fe/H] decreases. Nevertheless we note that the number of thick disk stars with a low error in age is significantly reduced in this set of ages. We note here that this temporal evolution in the thick disk would be more obvious if we were to consider the hαmr stars (or at least the oldest ones) as part of the thick disk. For example, in the works by Bensby et al. (2014); Haywood et al. (2013) the metallicity of the thick disk can reach values above solar.

Fig. 4 [Fe/H] as a function of ages from Gaia DR2 for the full sample (upper panel) and for the subsample of stars with errors in age lower than 1.5 Gyr. The different stellar populations are depicted with different colors and symbols as explained in the legend. The average and standard deviation of [Fe/H] in 1 Gyr size bins are shown with a black line and error bars.

5 Dependence of individual [X/Fe] ratios on [Fe/H] and age

In Fig. 7 we show the [X/Fe] ratios of the elements presented in DM17 as a function of age for the different populations in our sample together with a linear fit to thin disk stars. To complement the qualitative study of α elements wehave added the rederived oxygen abundances from Bertran de Lis et al. (2015), as shown in Fig. 4. A first look at this figure allows us to see the expected general trends. The ratios of α elements O, Mg, Si, Ca, and Ti respect to Fe show an increasing trend toward older ages with O and Mg showing the steepest trends (the slopes for thin disk stars are 0.026 and 0.02 dex Gyr⁻¹, respectively). These trends are in general similar to those observed in other works such as Nissen (2015), Bedell et al. (2018), Anders et al. (2018), Feuillet et al. (2018), Buder et al. (2019). Those elements are mostly produced in SNeII meanwhile Fe is produced mainly by SNeIa. Since the progenitors of SNeII are more massive than the progenitors of SNeIa, the ratios [α/Fe] will be higher at early ages in the Galaxy because massive stars have shorter lifetimes. However, the contribution of SNeII to O and Mg is higher than to the other α elements, which are also partially produced by SN Ia (e.g., Nomoto et al. 2013). Therefore, that could explain the steeper age trend for O and Mg. Since Al is mainly produced by core-collapse SNe (e.g., Andrews et al. 2017), the abundance-age trend is also similar to that of α elements, with a steep slope but a higher dispersion. However, we note that this positive trend might be limited by the lower limit in [Fe/H] of our sample because stars with [Fe/H] ≲ –1.5 dex (hence, old) show quite low [Al/Fe] values (see e.g., the compilation made by Prantzos et al. 2018). Finally, the element Zn also increases with age since it has an important contribution by neutrino winds during supernovae explosions of massive stars (e.g., Bisterzo et al. 2005, and references therein). Interestingly, this is the only element to show a similar behavior with age for both thin and thick stars and might represent a better global age proxy than [α/Fe] ratios.

It is also clear from Fig. 7 that thick disk stars present a stronger enrichment in α elements when compared to thin disk stars at similar age. This extra enrichment in thick disk stars is also observed for the r-process element Eu, which is mainly produced by neutron star mergers (Drout et al. 2017; Côté et al. 2018) and core-collapse supernovae (Travaglio et al. 1999). Both elements have massive progenitors, and in turn [Eu/Fe] also shows a rise toward older ages. Therefore, these trends support the results by Snaith et al. (2015) showing that the star formation rate was more intense in the thick disk than in the thin disk.

On the other hand, s-process elements show a decreasing trend of [X/Fe] as age increases. These elements are mainly produced in low-mass AGB stars so we can expect them to increase with time (for younger stars) due to the increasing and delayed contribution of low-mass stars as the Galaxy evolves. However, there seems to be a change of slope around 8 Gyr caused by thick disk stars, similar to the results by Battistini & Bensby (2016) for Sr and Zr⁷. The results of Spina et al. (2018) and Magrini et al. (2018) for Ba and Ce also seem to show a change of slope around 7 Gyr. We can observe that the abundance ratios of light-s elements (Sr, Y, Zr) are on average larger in the thick disk than in the thin disk for similar ages while the heavy-s elements (Ba, Ce) present similar enrichments at a given age. This means that the production of light-s elements was more efficient than that of heavy-s elements in the thick disk. Whether the overproduction of light-s elements at lower metallicities is produced by intermediate mass AGB stars or rotating massive stars is still unclear (see the discussion in DM17; Bisterzo et al. 2016; Prantzos et al. 2018).

In Fig. 7 we can also observe that heavy-s elements present a great dispersion at young ages whereas Sr and Y present tighter correlations with age. This seems to be caused by the different [X/Fe]-age relations depending on metallicity as explained in next section. Nd is also considered as a heavy-s element but it has a remarkable contribution of 44% from the r-process (Arlandini et al. 1999; Bisterzo et al. 2016), which is produced by massive stars (neutron star mergers or core-collapse supernovae), thus the behavior with age is balanced by the two sources and shows a rather flat trend. This is also the case of Cu, which has an important part produced through the weak-s process in massive stars during core He and shell C burning, where neutrons are provided by the ²²Ne(α, n)²⁵Mg reaction. Since the neutron source, ²²Ne, is originated from pre-existing CNO nuclei, it depends on the initial metallicity of the star. Therefore, the weak s-process is considered to be of secondary nature and that could explain why we do not see a steep rise toward older ages where stars are less metallic.

The relations of abundance ratios with [Fe/H] are represented in Fig. 8 with a color scale to show the ages. This plot manifestly represents the mixed influence of metallicity and age on GCE. For the α elements (Mg, Si, Ti), at a given [Fe/H] the abundances increase as age increases, as also shown in Fig. 5. The metal poor thin disk ([Fe/H] ≲ –0.3 dex) presents well mixed ages meanwhile the metal rich counterpart only has stars younger than ~5 Gyr. The hαmr stars are clearly older compared to thin disk stars at the same metallicity. The abundances of oxygen have much larger errors and thus the separation by ages is not as clear as for the other α elements. On the other hand, the Ca abundances do not show so clear separation between the different populations but still we can distinguish a general trend of increasing age and [Ca/Fe] as [Fe/H] diminishes.

A second group of elements is formed by Al, Cu and Zn which present a similar morphology. In this case the ages increase in a “diagonal” way, not only as the [X/Fe] ratio does but also with the drop of metallicity rather than at a fixed metallicity as happened for α elements. For the s-process elements Sr, Y, Zr and Ba we can see how the thin disk presents a very clear stratification of ages, also in a “diagonal” way, but with ages decreasing as both [X/Fe] and metallicity increases. However, the thick disk stars present a wider range of abundance ratios despite its narrow range of ages. The lower abundances of s-process elements in hαmr with respectto thin disk stars at the same [Fe/H] could be caused by their older ages. The peak of abundances of s-process elements around solar metallicity is formed by the youngest stars in the sample. However, the steep decrease of Ba abundances at supersolar metallicities does not depend on age. Ce and Nd present more mixed ages in the [X/Fe]–[Fe/H] plane, probably due to the larger errors of those abundance ratios. Finally, the abundances of Eu also show higher values as age increases and metallicity decreases, similar to O, but the ages are more mixed due to the larger uncertainties on the abundance derivation of this element.

Fig. 5 [α/Fe] as a function of [Fe/H] for the subsample of stars with errors in age from Gaia DR2 lower than 1.5 Gyr. The different stellar populations are depicted with different colors and symbols as explained in the legend.

Fig. 6 [α/Fe] as a function of [Fe/H] for the subsample of stars with errors in age from q2 lower than 2 Gyr. The different stellar populations are depicted with different colors and symbols as explained in the legend.

Fig. 7 [X/Fe] as a function of age for stars with an error in age smaller than 1.5 Gyr. The different stellar populations are depicted with different colors and symbols as explained in the legend. We note the different size of y axis for oxygen with respect to the rest of elements. The red line is a weighted linear fit to the thin disk stars to guide the eye on the general behavior of the trends. The coefs. values in each panel are the abscissa origin and the slope of the fit, respectively, together with the error (σ) of each coefficient.

Fig. 8 [X/Fe] as a function of [Fe/H] for stars with an error in age smaller than 1.5 Gyr. We note the different size of y axis for oxygen with respect to the rest of elements. The circles, triangles, squares and diamonds are the stars from the thin disk, thick disk, hαmr and halo.

6 Stellar age estimation using chemical abundances

As we have seen in previous sections stellar dating using chemical clocks is based on the different galactic chemical evolution of some species. The different contribution to the chemical evolution of the Galaxy of the SNeII, SNIa, and the low-mass AGB stars residuals opens the door to the stellar dating using certain surface chemical abundances (Nissen 2016). The work by da Silva et al. (2012) was the first exploring the relation with age of abundances ratios of Y or Sr over Mg, Al, or Zn. More recently, Nissen (2015, 2016) found that ratios of [Y/Mg], [Y/Al] or [Al/Mg] are precise age indicators in the case of solar twins stars. These are the so called chemical clocks and have been studied in other samples of solar twins (Spina et al. 2016; Tucci Maia et al. 2016) and very recently, in a bigger sample of stars within the AMBRE project (Titarenko et al. 2019). Moreover, these chemical clocks working over solar twin stars where confirmed using stars dated by asteroseismology (Nissen et al. 2017). However, Feltzing et al. (2017) and Delgado Mena et al. (2018) find that, when stars of different metallicities and/or effective temperatures are included, these simple correlations are not valid anymore. This can be clearly seen in Fig. 9 where the slopes of [X/Fe] ratios with age are shown for three metallicity groups. For instance, [Ca/Fe] has a mostly flat behavior with age for stars with –0.6 < [Fe/H] < –0.2 dex (black line) and –0.2 < [Fe/H] < 0.2 dex (blue line) but the metal rich stars ([Fe/H] > 0.2 dex, red line) present a negative trend. This is in contrast with the other α elements in most of the metallicity bins, which show an increase of abundance ratios with age. Recent studies report a new kind of supernovae subclass called calcium-rich gap transients (Perets et al. 2010) which can be an important contributor to the enrichment of Ca in the Galaxy, thus it might be the reason of the different evolution of Ca with time as compared to other α elements. Another interesting element is Cu, for which the age-abundance slope increases with metallicity. The case of most s-process elements is also remarkable, with a totally different slope for metal-rich and metal-poor subgroups, especially for Ce and Nd. The higher dispersion of heavy-s elements at younger ages as compared to light-s elements is certainly caused by the different metallicity of the stars. In other words, it seems that the nucleosynthesis channels producing heavy-s elements at younger ages have a stronger dependence on metallicity than those producing light-s elements. On the other hand, for several elements (O, Mg, Si, Ti, Zn, and Sr) the mix of stars with different metallicities can add dispersion to the abundance-age relation but without a significant change in slope.

To quantify the variations of the abundance-age slopes (for thin disk stars only) at different metallicity ranges we present in Fig. 10 the slopes of neutron-capture abundance ratios ([X/Fe]) with age, as a function of s-process contribution as shown for solar twins in Fig. 6 of Spina et al. (2018). In order to have enough stars with low errors in age within each range we have selected bins with a 0.4 dex width in [Fe/H] covering the full metallicity range of thin disk stars. The bin containing the solar twins (–0.2 < [Fe/H] < 0.2 dex) presents a similar result to that of Spina et al. (2018) except that Nd does not show a correlation with age for solar metallicity stars. These authors showed how the dependence with age, in other words the slope, is larger (and negative) as the s-process contribution increases, with r-process elements such as Eu presenting a negligible correlation with age. The authors conclude that the s-process production in the thin disk (where their solar twins belong) has been more active than the production of r-process elements. However, we note that this might not hold for all the metallicity ranges. At super-solar metallicities stars in the thin disk have produced less Ce and Nd in recent times and Y and Zr seem to be produced at a lower rate than for solar metallicity stars. On the other hand, Sr, presents a rather constant negative trend with age for all the metallicity bins. Nevertheless, we note that these results must be taken with caution because of the low significance of the slopes in such small subsamples.

Using the sample of reliable ages described in the previous sections, we have studied the most significant relations for estimating the stellar age. The main goals are: (1) find the dependences explaining the spread of [X/Fe] vs. age found when stars within a large range of metallicities and effective temperatures are observed; and (2) offer the best relations possible for a precise estimation of stellar ages. Among all the chemical species presented in DM17, there are two special cases: Eu (with large errors) and Zr I (only available for cool stars). There is a large number of stars with these abundances unknown so we have not used these two species to estimate stellar ages. In a similar way, we have not considered oxygen due to the large errors in the abundances and the difficulty of deriving it as compared to other α elements.

Fig. 9 [X/Fe] as a function of age for stars with an error in age smaller than 1.5 Gyr. We note the different size of y axis for oxygen withrespect to the rest of elements. The red, blue, and black lines are weighted linear fits to thin disk stars with [Fe/H] > 0.2 dex, –0.2 < [Fe/H] < 0.2 dex and –0.6 < [Fe/H] < –0.2 dex, respectively.

Fig. 10 Slopes (dex Gyr−1) of the [X/Fe]-age relations ([X/Fe] = a + b⋅Age) for neutron capture elements as a function of s-process contribution for thin disk stars att all metallicities and in four metallicity bins with a width of 0.4 dex.

6.1 Revising the chemical clocks candidates

The main characteristics of a chemical clock are its high correlation with age and, when a linear regression is fitted to it as a function of the stellar age, this regression is able to explain the main part of the observed variability. Most of the claimed chemical clocks in the literature fulfill these requirements. In this section we will analyze our best ages sampling trying to confirm and/or add new chemical clocks to those already known. To do so, we have studied the correlation between all the chemical species of our sample with age. In addition, and for reasons we describe in Sect. 6.2.2, we have also analyzed the correlation of the chemical species with stellar mass, T_eff and [Fe/H]. Thecorrelation has been determined using the Spearman correlation coefficient (ρ). It measures the rank correlation between two variables and so depicts monotonic relationships.

In Table 3 we show the values of ρ found. We will first focus on the correlations with the stellar age. To ensure the statistical significance of our results, we have also obtained the p-value of each determination. In all the cases, except for those correlations lower than 0.1, the p-value is < 0.001, that is, ρ has a high statistical significance. We find that [Zn/Fe], [Al/Fe], [Si/Fe], [Mg/Fe], [Ti I/Fe], [TiII/Fe], and [O/Fe] correlate with age (|ρ| > 0.7), [Ti II/Fe] a little better than [Ti I/Fe]. On the other hand, [Y/Fe] and [Sr/Fe] anticorrelate weakly with age, but they are the elements with the largest anticorrelation with age. Therefore, any ratio of [Y or Sr] over [Zn, Al, Si, Mg, or Ti] is a good candidate to be a chemical clock, that is, light-s elements over α elements (plus Zn and Al). This is something we could expect from the chemical elements formation history. Some of these chemical clocks have already been studied.

In Table 4 we show the correlation of those ratios and the stellar age. Some of the clocks are shown as a function of age in Fig 11. All chemical clocks candidates significantly anticorrelate with age with a p-value lower than 0.001. In other words, they are potentially good candidates for highly correlated 1D regressions with age, but some of them present a high dispersion caused by metallicity. We will verify this in the next section. With the information of Tables 3 and 4 we find some abundance ratios with a really high correlation with age, around |ρ| > 0.8. In the next sections we will study how to estimate stellar ages with abundance ratios. We also checked if we could get an improvement in the significance of the correlation with age⁸ by using a combination of α elements with Y abundances (e.g., [Y/(Mg+Si)]). The ρ values shown at the end of Table 4 are similar to those for chemical clocks with a single α element andfor simplicity we decided not to use them anymore.

For completeness, in Fig. 12 we show all the correlations found using the best ages subset (upper triangle). The value of the correlation represents the Spearman correlation coefficient (ρ). The color of the font represents the statistical significance of the value of ρ depending on the p-value, with red representing those with a p-value lower than 0.001, that is, the value has a very high statistical significance. All the correlations with |ρ| > 0.5 have a p-value <0.001. In this figure we also show the spread-plot of each pair of variables, with a fitting using a LOESS (local polynomial regression) curve to guide the eye (bottom triangle). As expected, the elements produced by similar processes/stars are highly correlated. For example, all the α elements and Eu are correlated among them with |ρ|≳ 0.7. Zn is highly correlated with Al, Mg, Si, and TiII since all those elements have an important contribution from SNeII. However, the s-process elements are not so correlated among them as the α elements, with |ρ| values closer to ~0.5. Finally, Nd is highly correlated with Eu, being both of them produced by the r-process. However, only ~50% of Nd is produced in such process, which takes place in massive progenitors. The other half is produced by the s-process in lower-mass stars. That would explain the fact that Nd is correlated with Ca, TiI and TiII (with |ρ| > 0.6) since those α elements are both produced in SNeII and SNIa (i.e., higher and lower mass progenitors).

Fig. 11 Some of the chemical clocks studied as a function of age with the metallicities in color scale. The weighted linear fits are the same as in Fig. 9.

6.2 Multivariable linear regressions for the estimation of stellar ages

One of the best options to condense all the information contained in a data sample in a simple relationship, with the aim of using it for stellar dating estimations in our case, is to obtain the most significant linear regressions possible. In fact, these relations try only to explain as much as possible the age variability. To avoid over-fitting and to keep the relations as simple as possible, we are going to present the study depending on the number of independent variables involved (or dimensions) up to the third dimension. We note that the linear regressions presented in next subsections are weighted by the errors on the variables (stellar parameters and abundance ratios).

Fig. 12 Correlations of all the variables with each other in the subset using the stars with the most accurate age determination (see text for details). In the upper triangular part, the Spearman statistic, ρ for each correlation is depicted. The color of the font represents the statistical significance of the correlations depending on the p-value, with red, black, and gray for values with p < 0.001, p < 0.01 and p < 0.1, respectively.In the lower triangular part we also show the spread-plot of each pair of variables, with LOESS fit just to guide the eye.

6.2.1 1D relations

These relations estimate the stellar age using only one independent variable: Age = f(X), where f(X) = (a ± Δa) + (b ± Δb) × X is the linear combination of the variable X. These are the relations we can usually find in the literature. In Table 5 we show the relations with an adjusted R² (adj-R²) > 0.65, where adj-R² is a measurement of the dependent variable variance explained by the independent variable X corrected by the number of dimensions involved. The p-value of all the adj-R² are, by far, lower than 0.001.

As expected, the relations with the largest adj-R² are those formed by the same abundance ratios with a large correlation with the age in Tables 3 and 4. [Y/Mg], [Sr/Mg], and [Y/Zn] are part of the chemical clocks proposed by da Silva et al. (2012), Nissen (2015, 2016), Nissen et al. (2017), Spina et al. (2016), Delgado Mena et al. (2018). On the other hand, Bensby et al. (2014), Haywood et al. (2013), Feltzing et al. (2017) speculated that Ti can also be a good species for a chemical clock but with certain scatter. Here we show that, in fact, [Y/TiII] and [Sr/TiII] are very good ones. We also evaluated the chemical clocks with the average of TiI and TiII (named as Ti in the tables). Since the stars included in the sample are not only solar twins, the adj-R² obtained are not impressive because of the behavior of abundances ratios at different metallicities (Feltzing et al. 2017; Delgado Mena et al. 2018), but values of adj-R² around 0.8 are remarkable.

To test the reliability of these regressions, we have randomly split the good ages subsample into a training and a control group, with a 70 and 30% of the stars respectively. We have re-obtained the regressions with the training group, estimated the ages for the testing group and compared them with the real ages (those in our catalog). In Table 5 we show the standard deviation (S.D.) of the difference estimated-real ages. For these best relations, the S.D. is on the order of 2 Gyrs. On the other hand, since an error of 1 Gyr (for example) is not the same for a 13 Gyr old star or for a 1 Gyr old star in relative terms, we have also calculated the mean relative difference between the estimated ages and the real ages. This value (represented as a fraction of one), is shown in the column “Mean rel”. In this case we can see that the mean error is, in fact, large for these 1D relations. The minimum error is on the order of 40%. A number of additional consistency tests for the regression models shown in this work can be found at the Appendix A. In Fig. 13 we also show a comparison example between the real ages and those obtained with the 1D formula using the [Y/TiII] chemical clock, only for a subset of stars with very low errors in age (<0.5 Gyr). The results, shown with yellow circles, demonstrate a clear dependence of the obtained ages with the stellar parameters [Fe/H] and T_eff, and the need to include more variables in the formulas to obtain the age (see next subsection).

We have also compared the relations substantiated by our database of solar twins with those found in the literature. To do so, we have selected a subset of stars with T_eff = [5677, 5877]K, and [Fe/H] = [−0.1, 0.1]. We allow the errors in age to be up to 2 Gyr to have enough stars for a meaningful comparison. In this case, the 1D relations are formulated as X = f(Age), so we have reproduced those in Nissen (2015, 2016), Nissen et al. (2017), Spina et al. (2016), Tucci Maia et al. (2016). Our results compared with those in the literature are shown in Table 6. All the coefficients in the literature and those we obtain are equivalent within uncertainties. We note, however, that we lack stars with low errors in age in the range 3–6 Gyr and this might explain the small differences in the regression coefficients.

Finally, our data sample support no correlation for the expression [Al/Mg] = f(Age) as shown in Nissen (2016). This result is reasonable since both [Al/Fe] and [Mg/Fe] correlate positively with age, and, therefore, its ratio is not expected to be a good chemical clock, as its the case when we rate a chemical element correlating with another anticorrelating with age.

6.2.2 2D relations

The next natural step was to analyze the improvement we obtain when we add a second independent variable to the relation. Taking into account the lessons in the literature (e.g., Feltzing et al. 2017) and in the previous sections of this work, one of the main sources for explaining the stellar age variability is the stellar metallicity. On the other hand, we cannot discard the impact of stellar structure and evolution on the transport of chemicals and hence on its chemical surface abundances (Salaris & Cassisi 2017; Dotter et al. 2017). For example, diffusion is considered to have a major effect on hotter stars in the turn-off (e.g., Bertelli Motta et al. 2018) Therefore, we have tested all the possible combinations Age = (a ± Δa) + (b ± Δb) * T_eff + (c ± Δc) * X, with X any of the chemical species or their ratios, and Age = (a ± Δa) + (b ± Δb) *[Fe/H] + (c ± Δc) * X.

Stellar mass is, in terms of the physics involved in the stellar chemical mixing with age, a good proxy for stellar evolution. Thus, we have also used M as an independent variable to obtain the age with a formula such as this: Age = (a ± Δa) + (b ± Δb) * M + (c ± Δc) * X. When two variables have a large correlation, their inclusion in a linear regression is not recommended since they provide redundant information.In Tables 3 and 4 we show the correlations of all our chemical elements with T_eff, [Fe/H], and M. We define that any pair of variables with a correlation with |ρ| > 0.7 should not be used simultaneously in the 2D formulas. In Table 7 we show the final relations obtained for all the combinations of chemical elements and T_eff with an adj-R² > 0.60. “Rel imp X or T_eff” represent the relative importance of the independent variables in the explanation of the observed variance, that is, the real impact of each variable in the estimation of the stellar age. In this table we can see:

• The inclusion of T_eff as additional variable increases adj-R² a 7.5%, going from 0.80 to 0.86 for the respective best cases. Therefore, with a chemical clock and the effective temperature as independent variables we can explain up to a 86% of the stellar age variance.

• The relative importance is almost always balanced, that is, both variables are needed for obtaining these results, with T_eff explaining around a 25% of the variance and the chemical clock the rest 75%.

• Seven of our proposed chemical clocks (all except [Sr/Al], [Sr/Zn], and [Y/Al]) have relations explaining at least a 64% of the stellar age variance.

• Some [X/Fe] ratios can also be a good chemical clock when used together with the effective temperature. This is the case of [Si/Fe], [Ti/Fe], and [TiII/Fe].

Therefore, we can confirm that stellar structure and evolution plays a role in the observed stellar age variance when compared with surface abundances. We would like to note the decreasing in the S.D., compared with Table 5. The 2D relations provide a S.D. which is around 10–20% more precise than 1D relations. Something similar happens in the case of the mean relative differences.

In Table 8 we show the best regression models obtained adding the stellar metallicity as second independent variable. In terms of the best statistics, the improvement reached compared with the 1D regressions is similar to that found in Table 7. The chemical clocks involved in the best relations are also similar as using T_eff as independent variables but the relative importance of [Fe/H] is higher than for T_eff. However, in this case, the best relations (ordered by adj-R²) are obtained with chemical clocks formed by Sr and Zn. We note that [Zn/Fe] and [Sr/Fe] already show a correlation with T_eff (see Table 3), so combining them with the [Fe/H] improves the determination of ages. On the other hand, the abundances ratios [Ti/Fe] and [Si/Fe], which have a strong correlation with [Fe/H] (see Table 3), are better combined with T_eff to obtain good 2D relations. Therefore, stellar metallicity also plays a role in the explanation of the stellar age variability when using chemical clocks.

Finally, using M as independent variable also provides similar results. The increase in terms of adj-R² with respectto the 1D relations is on the order of previous 2D relations and the S.D. values are lower than when using T_eff. The relative importances are even in a better balance than previously, probably because the stellar mass has a dependence on T_eff and [Fe/H]. In addition, these results show a possible link between stellar structure and surface abundances for explaining the observed variability of stellar ages. This points in the direction of different theoretical studies (Dotter et al. 2017; Salaris & Cassisi 2017) since stellar mass is a rough proxy for how surface abundances change with age.

6.2.3 3D relations

As a natural last step, we analyzed the improvements reached when taking into account [Fe/H], T_eff, and M in groups oftwo, together with a stellar clock for estimating stellar ages. We have used relations with the form Age = (a ± Δa) + (b ± Δb) * X + (c ± Δc) * Y + (d ± Δd) * Z, with X and Y two of [Fe/H], T_eff and M, and Z any of the abundance ratios. We have followed the same procedure as in the previous cases, with the best ages subset. In this case, the results show that for every Z, any combination of [Fe/H], T_eff and M provides almost the same results in terms of adj-R². As an example, in Table 10 we show the best results obtained in case of using [Fe/H] and T_eff as independent variables. Here we can see that the adj-R² is on the order (just slightly larger) of that of the 2D relations. The same results are found when using [Fe/H] and M, and T_eff and M. Our conclusion is that using a good proxy for describing stellar evolution, that counts for chemical mixing, contains almost all the information complementing that coming from the chemical clock. The remaining variability should come from another source and adding more general stellar characteristics does not add almost new information.

Besides these considerations, other conclusions of this table are:

• The best 3D relations present a remarkable adj-R² = 0.89.

• All the chemical clocks except [Sr/Al] and [Y/Al] (with an adj-R² close to our threshold) explain at least a 87% of the stellar age variance.

• An important consistency check of these expressions is that the relative importance of the different variables are in general a 40% for the chemical clock, a 20% for the effective temperature, and a 40% for [Fe/H], in other words, well balanced.

• The S.D. presents, in general, an improvement of a 30% with respect to the 1D relations. Something similar happens with the mean relative differences.

Therefore, as a general conclusion, a chemical clock in conjunction with one stellar parameter (either stellar mass or the effective temperature or the metallicity) or with the combination of both the effective temperature and the stellar metallicity can explain between 89 and 87% of the stellar age spread. This result confirms that stellar structure and evolution plays a role in stellar dating using surface chemical abundances, but we might need to consider additional physical processes, such as rotation, to understand the remaining variance. However, we note that most of our stars are slow rotators (v sin i ≲ 8 km s ⁻¹; Delgado Mena et al. 2015). Additionally, there are other factors not considered here, such as NLTE corrections or the atomic diffusion of elements with age that can affect the presented relations. Finally, we note that part of the spread observed in the [X/Fe] ratios is not astrophysical, and thus cannot be characterized. For example, the scatter in the abundance ratios has a dependence on the number of lines used for a given element (Adibekyan et al. 2015) and is, of course, affected by the errors in the stellar parameters. This is clearly shown by the low scatter of [X/Fe] values whenusing good quality spectra of solar twins with low errors in parameters (e.g., González Hernández et al. 2010).

6.2.4 Applicability

An ideal way to test if the formulas presented in previous subsections are reliable is to compare the ages we obtain from them with the real ages of the stars. However, our derived ages using isochrones can have large errors, especially in the case of cool stars. In Fig. 13 we show how the average differences between ages obtained with the formulas (for the [Y/TiII] clock) and those derived by isochrones depend on T_eff and [Fe/H]. In order to calculate the differences in ages we have used only those stars with errors in age lower than 0.5 Gyr, that is, stars for which the isochrone ages are the closest possible to real within our sample. From the upper panel in Fig. 13 it is clear that for metal-poor stars the differences can be quite large but for stars more metallic than –0.5 dex the ages retrieved from the empirical relation are quite close to the real age. The ages obtained with the 3D formula (blue circles) are in better agreement with the derived ages in most of the metallicity bins. On the other hand, if we look at the applicability as a function of T_eff (bottom panel), it is obvious that the 3D formula (which uses T_eff) does not work well for cool stars. This is due to the fact that the linear regressions have been calculated using stars mostly hotter than 5300 K. Therefore, we do not recommend using those formulas for cool stars, nor for stars hotter than 6500 K due to the low number of hot stars used to derived the formulas. In a similar way, we do not recommend to use the formulas for metal poor stars since the regressions are obtained using stars mostly with [Fe/H] > –0.8 dex. In most of the T_eff bins the use of 3D clocks provides the most similar ages except for the stars around 6000 K for which the formulas tend to overerestimate the age.

In Fig. 14 we have compared the results obtained by using 3D formulas with different abundance ratios (and T_eff and [Fe/H]). The differences and dispersions in a given bin are very similar among the different formulas, hence the different formulas of the same dimension provide comparable results. However, when deriving ages for a given star using formulas with different dimensions or parameters the results might not be similar. As a rule, we would recommend use of the 3D formulas or the 2D formula with the stellar mass, since they provide the best results. From the bottom panel of Fig. 14 we see that the ages from those three formulas are slightly better than the age from formula with [Y/TiII] for stars around 6000 K (bottom panel of Fig. 13). Therefore, if several abundances ratios are available one can choose to get the ages from different formulas and get an average value. We also recall that when using these formulas for a given sample, the errors on the abundances of such sample (together with the errors in the coefficients) must be considered in order to have realistic errors on the determined age.

Fig. 13 Comparison of the [Y/TiII] 1D, 2D (using [Fe/H]) and 3D clocks with the Padova ages for a small sample of stars with error in age lower than 0.5 Gyr. Average differences between age from formula and real age in different [Fe/H] bins (upper panel) and different Teff bins (lower panel). The dispersion of the differences are shown with the error bars when there are more than one star in a given bin.

Fig. 14 Comparison of the 3D formulas using the chemical clocks [Y/Zn], [YMg], and [Sr/Ti] with the Padova ages for a small sample of stars with error in age lower than 0.5 Gyr. Average differences between age from formula and real age in different [Fe/H] bins (upper panel) and different Teff bins (lower panel). The dispersion of the differences are shown with the error bars when there are more than one star in a given bin.

7 Summary and conclusions

In this workwe derived different set of ages using either parallaxes from Gaia DR2 and HIPPARCOS or spectroscopic log g for a sample of more than 1000 solar neighborhood stars belonging to the HARPS-GTO program. The chemical abundances of those stars were presented in previous works. The aim of this work is twofold. On the one hand we evaluated how abundance ratios of elements with different nucleosynthetic origin evolve with time. On the other hand we have provided different empirical relations to determine ages from abundance ratios and stellar parameters. The main results of this work can be summarized as follows.

• Our results confirm the large dispersion of the age-metallicity relation in the solar neighborhood, which increases with age except for stars older than 12 Gyr, belonging to the thick disk. As found in previous works, the most metallic stars in our sample are not young, indicating that they had time to migrate to the solar neighborhood. The use of [α/Fe] (being α the average of Mg, Si, and Ti) or [O/Fe] provides a much tighter relation with age which becomes steeper for thick disk stars. On the other hand, we found that [Zn/Fe] presents a single tight relation with age valid for both thin disk and thick stars. This is because thick disk stars are less enhanced in Zn than α elements with respect to thin disk stars.

• We find that the ages of thin disk stars in our sample reach up to 11 Gyr in few cases. The oldest thin disk stars have a similar age and [α/Fe] content as the oldest hαmr but with a lower metallicity. These old thin disk stars have been regarded as coming from the outer disk (Haywood et al. 2013). Stars of a given [Fe/H] increase their ages as they increase their [α/Fe], thus, hαmr are well differentiated from thin disk stars of similar metallicities. When looking at s-process elements, the ages decrease with both [Fe/H] and [X/Fe] increasing. The maximum peak of [X/Fe] ratios for s-process elements is formed by the youngest stars which are those having solar metallicity.

• Thick disk stars present a stronger enrichment in α elements when compared to thin disk stars of similar age. This is also true for the r-process element Eu and for the light s-process elements. However, heavy s-process elements show a lower level of enrichment in the thick disk compared to light-s elements. Intermediate mass AGB stars or rotating massive stars have been identified as possible responsibles for the light-s elements overproduction in the thick disk.

• The abundances of α elements, Al, Zn, and Eu with respect to Fe generally increase with age whereas s-process elements over Fe diminish with age. The elements with a contribution of both massive stars and lower mass stars (such as Cu and Nd) show a rather flat behavior with age. However, the trends present large dispersions for some of the abundance ratios, mostly caused by the wide range of metallicities in each age bin. In addition, the slopes of some of the [X/Fe]-age relations change with [Fe/H] indicating the strong influence of metallicity for some nucleosynthesis channels. We find that such slopes have a remarkable change for metal-rich stars in the cases of [Ca/Fe], [Cu/Fe], [Ce/Fe] and [Nd/Fe]. The slopes of [Y/Fe], [ZrII/Fe] and [Ba/Fe] vs age become rather flat for metal rich stars whereas some ratios such as [Mg/Fe], [Si/Fe], [Ti/Fe], [Zn/Fe] and [Sr/Fe] have a quite constant trend with age regardless of the metallicity.

• The observed variation of abundance-age slopes proves that the use of simple linear functions to derive ages from certain abundance ratios (also called chemical clocks) is limited to certain ranges of metallicity. Therefore, we investigate how the inclusion of stellar parameters in multivariable linear regressions can help in estimating stellar ages. Our results shows that by using different chemical clocks combined with one or two stellar parameters (T_eff, [Fe/H] or stellar mass) we can explain up to a 89% of the age variability. The derived formulas can thus be used as an age-proxy for stars for which the derivation of stellar ages through other methods is not possible. We note, however, that the empirical relations presented in this work have limitations and should not be used for stars outside the parameters range of our sample.

The overall results of this work show how important is to add the stellar age information when studying the GCE although we must be aware of the uncertainties involved in the derivation of such ages. This is especially a disadvantage for cool stars, because the large errors on stellar parameters will always prevent to get precise ages when using the isochrone method. Future asteroseismic observations will provide more reliable ages for different kind of stars. This will allow testing and revising of the results presented here and recalibration of the formulas when necessary, making them useful in a broader parameter range.

Acknowledgements

E.D.M., V.Zh.A., M.T., N.C.S., and S.G.S. acknowledge the support from Fundação para a Ciência e a Tecnologia (FCT) through national funds and from FEDER through COMPETE2020 by the following grants UID/FIS/04434/2013 & POCI–01–0145-FEDER–007672, PTDC/FIS-AST/7073/2014 & POCI–01–0145-FEDER–016880 and PTDC/FIS-AST/1526/2014 & POCI–01–0145-FEDER–016886. E.D.M. acknowledges the support by the fellowship SFRH/BPD/76606/2011 funded by FCT (Portugal) and by the Investigador FCT contract IF/00849/2015/CP1273/CT0003 and in the form of an exploratory project with the same reference. V.Zh.A., N.C.S. and S.G.S. also acknowledge the support from FCT through Investigador FCT contracts IF/00650/2015/CP1273/CT0001, IF/00169/2012/CP0150/CT0002 and IF/00028/2014/CP1215/CT0002 funded by FCT (Portugal) and POPH/FSE (EC). A.M. acknowledges funding from the European Union’s Horizon 2020 research and innovation program under the Marie Sklodowska–Curie grant agreement No 749962 (project THOT). W.J.C. acknowledges support from the UK Science and Technology Facilities Council (STFC). Funding for the Stellar Astrophysics Centre is provided by The Danish National Research Foundation (Grant agreement no.: DNRF106). J.I.G.H. acknowledges financial support from the Spanish Ministry project MINECO AYA2017-86389-P, and from the Spanish MINECO under the 2013 Ramón y Cajal program MINECO RYC-2013-14875. A.C.S.F. is supported by grant 234989/2014-9 from CNPq (Brazil). Finally, we thank the anonymous referee, whose detailed review helped to improve this paper. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. This research has made use of the SIMBAD database operated at CDS, Strasbourg (France). The core analysis was performed with R version 3.3.1 (R Core Team 2016), RStudio Version 1.0.143, and the R libraries dplyr 0.5.0 (Wickham & Francois 2016), foreach 1.4.3 (Analytics & Weston 2015), corrplot 0.86 (Wei & Simko 2017), parallel (R Core Team 2016), RcppParallel (Allaire et al. 2016), and relaimpo 2.2-2 (Grömping 2006).

Appendix A Consistency of the regressions

We have developed a large number of tests to ensure the consistency and reliability of the regressions shown at Sect. 6.2. We have selected randomly 36 of the relations shown in Tables 5–10 for illustrating these tests.

In Fig. A.1 we show quantile–quantile plots of all the selected relations in a form to show whether the standardized residuals are normally distributed. The ordered standardized residuals are plotted on the ordinate of each plot, while the expected order statistics from a standard normal distribution are on the abscissa. Points close to the straight line are consistent with a normal distribution. The different classification of the stars as being part of the thin disk, thick disk, halo, or hαmr is shown in different colors. In this figure we can see that all the relations, in general, follow this straight line. Therefore, the use of linear regression is justified, but we must explore the special case of relations taking into account only thin disk stars.

In Fig. A.2 we present the residuals as a function of the fitted values of these relations. Any clear trend in these residuals can be a signature of inaccurate or inefficient regression. We have also added a LOESS (local polynomial regression) curve to guide the eye. In this figure we can see that, in general, there are no clear trends in the distributions of the residuals. In every plot, the main body of points is randomly distributed around the value zero. We can only seeboundary effects, where for the extreme age cases (close to zero and close to the limit of 13 Gyrs) the departure from the regressed model are larger. In addition, we can see in some cases a number of boundary stars with a large residual.

Attending to the histogram of these residuals (Fig. A.3), we find that the distributions are mostly close to Gaussian-like distributions, fulfilling one of the main assumptions of regressions models.

Fig. A.1 Quantile–quantile plots of the relations. See text for details.

Fig. A.2 Residuals vs. fitted values of the relations. The line is a LOESS curve to guide the eye. See text for details.

Fig. A.3 Histogram of residuals of the relations. See text for details.

References

All Tables

All Figures

	Fig. 1 Normalized distribution (per stellar population) of ages for the full sample (upper panel) and a subset with small errors in age (lower panel) using Gaia DR2 parallaxes. The yellow, red, blue and purple histograms represent the thin, hαmr, thick and halo stars, respectively.
In the text

	Fig. 2 Comparison of ages and masses obtained using Gaia DR2 or HIPPARCOS parallaxes and obtained with package q2. The red dots are the stars with errors in age lower than 1.5 Gyr. The green dashed lines shows the zero differences and the blue lines show the mean and standard deviation.
In the text

	Fig. 3 Comparison of age differences vs. parallax differences between Gaia DR2 and HIPPARCOS. The V magnitudes are shown in a color scale.
In the text

	Fig. 4 [Fe/H] as a function of ages from Gaia DR2 for the full sample (upper panel) and for the subsample of stars with errors in age lower than 1.5 Gyr. The different stellar populations are depicted with different colors and symbols as explained in the legend. The average and standard deviation of [Fe/H] in 1 Gyr size bins are shown with a black line and error bars.
In the text

	Fig. 5 [α/Fe] as a function of [Fe/H] for the subsample of stars with errors in age from Gaia DR2 lower than 1.5 Gyr. The different stellar populations are depicted with different colors and symbols as explained in the legend.
In the text

	Fig. 6 [α/Fe] as a function of [Fe/H] for the subsample of stars with errors in age from q2 lower than 2 Gyr. The different stellar populations are depicted with different colors and symbols as explained in the legend.
In the text

Fig. 7 [X/Fe] as a function of age for stars with an error in age smaller than 1.5 Gyr. The different stellar populations are depicted with different colors and symbols as explained in the legend. We note the different size of y axis for oxygen with respect to the rest of elements. The red line is a weighted linear fit to the thin disk stars to guide the eye on the general behavior of the trends. The coefs. values in each panel are the abscissa origin and the slope of the fit, respectively, together with the error (σ) of each coefficient.

In the text

	Fig. 8 [X/Fe] as a function of [Fe/H] for stars with an error in age smaller than 1.5 Gyr. We note the different size of y axis for oxygen with respect to the rest of elements. The circles, triangles, squares and diamonds are the stars from the thin disk, thick disk, hαmr and halo.
In the text

	Fig. 9 [X/Fe] as a function of age for stars with an error in age smaller than 1.5 Gyr. We note the different size of y axis for oxygen withrespect to the rest of elements. The red, blue, and black lines are weighted linear fits to thin disk stars with [Fe/H] > 0.2 dex, –0.2 < [Fe/H] < 0.2 dex and –0.6 < [Fe/H] < –0.2 dex, respectively.
In the text

	Fig. 10 Slopes (dex Gyr−1) of the [X/Fe]-age relations ([X/Fe] = a + b⋅Age) for neutron capture elements as a function of s-process contribution for thin disk stars att all metallicities and in four metallicity bins with a width of 0.4 dex.
In the text

	Fig. 11 Some of the chemical clocks studied as a function of age with the metallicities in color scale. The weighted linear fits are the same as in Fig. 9.
In the text

Fig. 12 Correlations of all the variables with each other in the subset using the stars with the most accurate age determination (see text for details). In the upper triangular part, the Spearman statistic, ρ for each correlation is depicted. The color of the font represents the statistical significance of the correlations depending on the p-value, with red, black, and gray for values with p < 0.001, p < 0.01 and p < 0.1, respectively.In the lower triangular part we also show the spread-plot of each pair of variables, with LOESS fit just to guide the eye.

In the text

Fig. 13 Comparison of the [Y/TiII] 1D, 2D (using [Fe/H]) and 3D clocks with the Padova ages for a small sample of stars with error in age lower than 0.5 Gyr. Average differences between age from formula and real age in different [Fe/H] bins (upper panel) and different Teff bins (lower panel). The dispersion of the differences are shown with the error bars when there are more than one star in a given bin.

In the text

Fig. 14 Comparison of the 3D formulas using the chemical clocks [Y/Zn], [YMg], and [Sr/Ti] with the Padova ages for a small sample of stars with error in age lower than 0.5 Gyr. Average differences between age from formula and real age in different [Fe/H] bins (upper panel) and different Teff bins (lower panel). The dispersion of the differences are shown with the error bars when there are more than one star in a given bin.

In the text

	Fig. A.1 Quantile–quantile plots of the relations. See text for details.
In the text

	Fig. A.2 Residuals vs. fitted values of the relations. The line is a LOESS curve to guide the eye. See text for details.
In the text

	Fig. A.3 Histogram of residuals of the relations. See text for details.
In the text