# Double-Dipping: A New Relation between Stellar Rotation and Starspot
Activity

Gibor Basri and Hieu T. Nguyen  Astronomy Department, University of
California, Berkeley, CA 94720<EMAIL_ADDRESS>
###### Abstract

We report the discovery of a new relationship between a simple morphological
characteristic of light curves produced by starspots and stellar rotation
periods. The characteristic we examine is whether the light curve exhibits one
dip or two during a single rotation. We analyze thousands of Kepler light
curves of main sequence stars from 3200-6200K. Almost all the stars exhibit
segments of their light curve that contain either single or double dip
segments (very few have more than two significant dips per rotation). We
define a variable, the “single/double ratio” (SDR) that expresses the ratio of
the time spent in single mode to the time spent in double mode. Unexpectedly,
there is a strong relationship between the SDR and the stellar rotation
period, in the sense that longer periods come with a larger fraction of double
segments. Even more unexpectedly, the slopes of the SDR-Period relations are a
clear function of stellar temperature. We also show that the relationships of
spot variability amplitude ($R_{var}$) to rotation period have similar levels
of scatter, slopes, and dependence on temperature as the SDR-Period relations.
Finally, the median $R_{var}$ of single segments tends to be about twice that
of double segments in a given light curve. We offer some tentative
interpretations of these new results in terms of starspot coverage and
lifetimes. It will be fruitful to look further into this novel “rotation-
activity” relation, and better understand what information these aspects of
light curve morphology bring to our knowledge of stellar magnetic activity.

starspots — stars: magnetic field — stars: activity — stars: late-type —
stars: solar-type — stars: rotation

††slugcomment: to appear in The Astronomical Journal

##  1 Introduction

Sunspots have fascinated scientists since being documented by Galileo in the
early 17th century. There have been centuries of observations of their
behavior on the Sun, yielding insights into its surface and internal rotation
and the operation of its magnetic dynamo. Photometric changes in other stars
ascribed to starspots have also been observed for more than a century, but
with vastly poorer time coverage and precision and little information on their
sizes and positions. Things got slightly better with the advent of Doppler
Imaging (Vogt & Penrod, 1983) and more recently Zeeman Doppler Imaging (Semel,
1989), which began to provide information on the size and position of spots
(e.g. Strassmeier (2009); See et al. (2017)).

With the advent of precision space photometry, particularly from the COROT
(Baglin, 2003) and Kepler (Borucki et al., 2010) space telescopes, this
situation has been markedly improved. In particular, the Kepler mission
provides nearly continuous coverage over 4 years for well over a hundred
thousand stars, with sufficient precision to detect individual sunspots if the
Sun were hundreds of parsecs away. This has yielded a wealth of information on
stellar rotation periods (e.g. McQuillan, Mazeh & Aigrain (2014)), and has
been mined to learn more about spot coverage and variability, and potentially
differential rotation for thousands of stars. There have been many other
detailed analyses of individual stars with various types of mapping techniques
employed, but here we are concerned with statistical behavior of groups of
stars. The work that is closest in spirit to this paper is an analysis of
about a thousand Kepler stars by Arkhypov et al. (2018). Although their sample
is much smaller and restricted to the lower half of our period distribution,
they also consider the behavior of the stars per rotation period, and make
some use of period harmonics. Another paper with a restricted sample size (but
only 2 discrete periods) is that of Giles, Collier Cameron & Haywood (2017).
They concentrate on the question of spot lifetimes, and utilize
characteristics of the entire light curve, but are also affected by some of
the issues of concern here. This paper takes a much simpler approach (which is
much more easily visualized) and considers a much larger sample with the full
period range. Some important issues remain for all work to date on how to
translate Kepler pipeline differential intensities into correct understandings
of starspot distributions and evolution, and the actual information content of
broadband light curves (particularly a better understanding of the
degeneracies they hide).

Starspots are, of course, one manifestation of stellar magnetic activity. One
of the most basic facts about this activity is that it depends on stellar
mass, age, and rotation period. A variety of “rotation-activity” relations
have been studied over the past several decades (e.g. Noyes et al. (1984);
Reiners, Schüssler & Passegger (2014)). They have concentrated on diagnostics
of magnetic heating, from the chromosphere to the corona. In general, the
faster a star with a convective envelope and a given mass is rotating, the
stronger is the magnetic activity and the brighter its emission diagnostics.
It is also true that this activity induces magnetic braking in most cases,
which causes stars to spin down over time and thus their activity to decrease.
This leads to the method of gyrochronology, which relates age, mass, and
rotation (e.g. Meibom et al. (2015)). Much less work has been done on how
starspots fit into this picture because of the lack of appropriate data until
recently. General relations between the amount of starspot variability and
stellar rotation have been shown, for example, by Reinhold et al. (2013);
McQuillan, Mazeh & Aigrain (2014); Giles, Collier Cameron & Haywood (2017);
Arkhypov et al. (2018). This paper looks into a barely explored aspect of
starspot signatures, and develops a new way of gaining information based on
the changing shapes of stellar light curves over time.

###  1.1 Kepler light curve morphologies

In examining Kepler light curves, one of the most obvious features is that a
given light curve often contains segments with systematic changes of amplitude
over time. One can determine the rotation period of a star if it has an
asymmetric spot distribution that retains some features for several rotations
and produces a measurable signal. It is then possible to characterize
properties of the light curve over each rotation. One of the simplest of these
is a count of how many dips occur during a rotation. It turns out that because
of the poor spatial resolution inherent in the light curve by itself, there
are typically only one or two dips per rotation. It is important to keep in
mind that although it is tempting to associate a dip in the light curve with a
spot of a given size, the reality is that it reflects the sum of spots over
most of a hemisphere. Thus, when we refer to a “spot”, we are generally
talking about a number of spots that have a spatial correlation with each
other (which might be strong or fairly vague). Nonetheless, we will continue
to use the term “spot” for convenience.

For two long-lived spots of equal size rotating with different periods, the
light curve will have a larger amplitude when they are relatively near each
other on the stellar surface (generating a single dip while in view during a
rotation) but a smaller amplitude with double dips per rotation when they are
on opposite sides of the star and not reinforcing each other. The advantage of
integrating over a whole rotation period is that (unless spots grow and decay
on the timescale of a rotation) all potentially visible spots will produce a
signal for part of the rotation, regardless of how they may be distributed.

Another possibility is that the changes in amplitude of the light curve are
due to spot evolution, with larger light deficits being produced when there
are more or larger spots present. One would have to also invoke preferentially
active longitude bands or long spot lifetimes (and no differential rotation)
if these spots are to preferentially produce single or double dips over many
rotations, since the hemispheric distribution would have to be maintained at
some level. It is less obvious in this case what the relation between single
and double dip amplitudes should be. In practice, it will often be the case
that both differential rotation and spot evolution are important. For this
paper it does not matter what is producing the changes in the light curve; we
are examining the morphologies themselves.

##  2 A New Diagnostic: the Single/Double Ratio

One of the most noticeable characteristics of Kepler light curves of main
sequence stars is that many of them show obvious periodicity (when they show
measurable variability at all). This may not seem surprising until one
remembers that the Sun itself does not usually present an obviously periodic
light curve. The methods used to find stellar rotation periods on Kepler stars
fail in significant chunks of the solar light curve (despite its excellent
signal-to-noise). However a number of authors have been able to determine
rotation periods for tens of thousands of stars in the Kepler sample. An oft-
quoted collection of them is in McQuillan, Mazeh & Aigrain (2014), and we
utilize that paper both as a source of targets and for adopted rotation
periods. It is worth pointing out here that there is an even larger body of
main sequence stars which do not display a measurable rotation signal. These
are likely to be older stars with too little activity, so the conclusions of
this paper really only apply to stars that are above that threshold of
activity.

McQuillan, Mazeh & Aigrain (2014) used an auto-correlation method, and noted
the common appearance of two repeating but dissimilar correlation function
features which occur with a certain time spacing, and about half that spacing.
This leads to the question as to which spacing represents the true stellar
rotation period (the presence of two periods can also confuse periodogram
methods of finding periods). McQuillan, Mazeh & Aigrain (2014) noted that
generally one of the features looks weaker than the other, and that pattern
often repeats for a while. They used this asymmetry to distinguish between the
case of two features per period or just one, allowing an assignment of a true
period (rather than its first harmonic).

Obviously, it is hard to claim there are two dips per rotation unless there is
some sort of repeating asymmetry between them; if they all looked the same or
fluctuated randomly then the assigned period would likely be the spacing
between adjacent dips. Another good way to distinguish single from double dips
is to look for the largest repeating separation between dips (widest light
curve features). That is because it is not really plausible for a star to
require two rotations to produce what looks like a single dip; the star would
have to get steadily darker over one rotation then steadily brighter in a
similar fashion over the next rotation. This is especially unlikely if the
single dips repeat several times. These points suggest further techniques for
distinguishing between true and half periods by keeping track of the
separations between dips as well as the amplitudes of successive dips.

The basic morphological question we ask is: what percentage of the time does a
star display a single dip per rotation, and what percentage of the time are
there double dips per rotation? It is this characteristic of the light curves
that is the primary point of this paper. We specifically define the quantity
of interest as the logarithm of the ratio of the time spent in single dip mode
to the time spent in double dip mode (which we will refer to as the
single/double ratio, or SDR). This ratio is evaluated over the whole Kepler
record for each star, typically a four year timespan. Thus an SDR of 0.0 means
that the star spent as long in the single dip mode as in the double dip mode,
an SDR of -1.0 means that it spent ten time longer in the double mode, and an
SDR of 1.0 means it spent ten times longer in the single mode. The main
results in this paper are that the SDR is a strong function of the rotation
period of the star, and the slope of the relation between SDR and rotation
period is a clear function of the effective temperature of the star.

Figure 1: Portions of a particular star’s light curve that contains both
single and double dip segments. The top panel shows a sample single dip
segment (solid) and double dip segment (dashed) on the same relative
timescale. The actual starting days are 1030 and 480 respectively. The
ordinate is the differential intensity in parts per thousand (ppt). The bottom
panel’s upper plot shows a longer section of the full light curve. For that
the ordinate is also in differential ppt, but with zero displaced up to the
solid line. The asterisks below mark the separation between successive dips;
for those the ordinate is in days (with zero at the bottom). The dashed line
indicates 75% of the rotation period of 9.94 days, the boundary below which we
classify dips as double. These lower points straddle half the period, and
sometimes display coherent temporal structure.

The data used in this analysis are PDC-SAP pipeline products from 2015,
obtained from the Kepler Archive at the MAST. We checked a number of cases
from the last release of this pipeline product (DR25) and found it is not
necessary for the purposes of this paper to re-download the whole dataset. The
long cadence Kepler data have a 29 minute time resolution, which is far finer
than needed to track the variability of spots. To speed things up and to
eliminate features that are irrelevant to our purpose (because they are too
fast to be produced by spots), we first rebinned the data by a factor of 10 to
timesteps of about 0.2 days. To further reduce the effect of “fake dips” due
to noise instead of spots, we implemented boxcar smoothing in which the
smoothing width for a given light curve is one eighth of the period found by
McQuillan, Mazeh & Aigrain (2014). This choice is heuristic; it means that the
dips being counted are visibly significant and have timescales that are
plausibly due to starspots. Our results are not qualitatively sensitive to the
particular choice of smoothing.

We locate the local maxima(peaks)/minima(dips) by comparing the relative flux
at any given point in time with its neighbors, point by point throughout the
whole light curve. We required that an extremum be present in a comparison of
the nearest 4 points (to further eliminate too fast features). It turns out
not to matter much whether one tracks dips or peaks, but we decided to use
dips as our main focus since they are directly related to starspot presence.
The question of small dips that are real but not counted here, and inflections
that are not quite dips, will be taken up in a subsequent paper that looks at
the behavior of these features in amplitude and phase over time (to address in
detail the question of whether they are likely indicators of spot drifting).

After finding the time at which each dip occurs, we calculate all the temporal
separations from one dip to the next. Then we compare these separations to the
rotation period, to sort out the single dips and double dips. If the
separations are greater than 75% of the rotation period, the segments of light
curve containing them are defined to be single dip segments; otherwise they
are considered double dip segments. Figure 1 shows an example of what we are
measuring. This star switches between single and double dip modes, often
spending a few rotations in each (though not always). The single separations
are mostly near the rotation period of 9.94 days; the double separations are
more scattered around half that period. Structures that might be associated
with phase drift are seen in the behavior of the some of the double dip
separations (near 770 or 910). These tend to switch high and low as a
secondary dip is closer to the previous and further from the next primary dip.

The choice of the boundary timescale midway between the rotation period and
half the period is motivated in part by looking at histograms of the
separations in light curves with significant segments of both single and
double dips. One tends to see a concentration of separations around the period
and half period, with a minimum about midway between them (Figure 2).
Sometimes the separations approach our (somewhat arbitrary) boundary, but
those are usually clearly associated with double dips. Of course, the
histograms look different from star to star, and sometimes there is not so
clear a separation between major histogram peaks (or even 2 major peaks at
all). We are thus forcing a binary characterization of more complicated
distributions, which should be kept in mind. However as we describe below, the
stars behave in a very systematic way despite this simplification.

Finally, we add up the durations of all the segments that have been
characterized as either single or double to find the total duration that the
star spends in single/double mode, and thereby determine the SDR. In a future
paper we will analyze the behavior of the timing of dips (along with their
amplitudes) in more detail, as this appears to be one fruitful approach to
understanding the differences between differential rotation and evolution.

Figure 2: A few example histograms of dip separations. The top row has stars
that are predominantly double dip, the middle row has about equal durations of
single and double mode, and the bottom row has predominantly single dips. In
each case the vertical line indicates 75% of the rotation period; the point of
division between single and double dips. Each case shows the minimum in the
number of dip separations that occurs near this chosen boundary. Each title
has the Kepler ID and the rotation period in days.

##  3 Analysis and Discussion

McQuillan, Mazeh & Aigrain (2014) provide a well-determined and widely cited
compilation of rotation periods for slightly over 34,000 stars from the Kepler
target list. Most of these are main sequence stars; we restricted our sample
by demanding that the log of the gravity obtained from the 2016 compilation of
stellar parameters in the Exoplanet Archive (Akeson et al., 2013) was greater
than 4.3. This was also the source of our stellar temperatures. To avoid stars
whose periods are due to pulsations rather than starspots, we restricted our
period range to be greater than 3 days, and the temperatures to be between
3200-6200K. This produced an overall sample of 26,628 stars, which we then
broke into 15 temperature groups of 200K each (Table 1). The groups do not
have the same number of stars in them due to the way that the Kepler team
chose targets; the coolest groups are particularly small and may not be fully
statistically representative of their population. The period ranges for the
cooler star samples are larger than for warmer stars because those samples
contain more long period cases (e.g. McQuillan, Mazeh & Aigrain (2014)). In
addition, we cut out the longest 2.5% of periods in each group, because those
longest period tails have much smaller densities of stars. Those tails all
have flat SDRs with period (in the double regime), and the rotation periods
are subject to more uncertainty because of potential problems with the
pipeline reductions for such long timespans.

Figures 3 through 5 display the main results of this paper. Each panel shows
the relationship of SDR to rotation period for one of the temperature groups;
we do not show the smaller groups below 3800K. Perhaps the first thing to note
is that the points lie almost entirely between -1.0 and 1.0, which means that
each of the light curves spends at least a tenth of the time in single dip
mode and similarly in double dip mode. The points which have an SDR of 0.0
spend equal time in each state. Such cases are best for studying the ways in
which light curves transition from single to double and back in detail
(although we leave that for a later paper). It is also of note that the
majority (between two-thirds and three-quarters in each group) of the stars
have a negative SDR, which means that the double dip mode is preferred by most
stars (column 3 of Table 1).

Figure 3: The relation between the SDR and rotation period in days for the
warmest temperature groups. Note the changing period ranges in the different
panels. A cubic fit to the sample is shown in red (note that SDR is a
logarithmic quantity). A linear fit to each group of stars is also shown in
green. Note that except for the warmest stars, the cubic fit tends to flatten
a bit at either end. The purple points were not used for the final fits (they
were outliers in initial cubic fits). The quiet and active Sun is indicated
with the yellow diamonds.  Figure 4: The same as Fig. 3 but for the middle
temperature groups. The cubic fits tend to flatten a bit at either end.
Figure 5: The same as Fig. 3 but for the coolest temperature groups. The cubic
fits for the two coolest groups rise at the short period end (unlike warmer
groups).

Each panel shows a very obvious (inverse) correlation between SDR and rotation
period, with the faster rotators tending to be more single and the slower
rotators tending to be more double. We determined a simple linear fit between
SDR (which is a logarithmic quantity) and period, and in each case it provides
a good description of the trends. It is very interesting that both the slopes
and the fraction of the points with positive SDR have systematic behaviors
with stellar temperature, as shown in Fig. 6. The slope steepens from
3200-5800K then starts to become shallower again for hotter stars. The steeper
slopes also tend to go with a greater proportion of single dip light curves
(except at the cool end). The values of the slopes and fractions are tabulated
in Table 1.

We also computed cubic fits to the points to capture potential curvature in
the relations. This was done iteratively: after the first fit we discarded
points with residuals of more than 2.5$\sigma$ before fitting again. There is
a tendency for the middle temperature groups to show flattening at either end.
The cubic fits provide comparability with the variability relations discussed
next, which showed more obvious curvatures. It is clear for the SDR-Period
relations, however, that linear fits are quite adequate.

Figure 6: The systematic behavior of the slopes of the linear fits in Figs. 3
through 5 with stellar temperature. Note that the negative of the slopes is
shown; higher numbers indicate steeper slopes. Also shown is a representation
of the fraction of the stars that have positive SDR, which is also systematic
with temperature. For comparison, the slopes of the long period sides of the
$R_{var}$-Period relations (discussed below) are also shown. They are quite
similar to the slopes of the SDR-Period relations for most temperatures. We
have adjusted the scales for each variable so they can be easily compared by
eye in this Figure.

###  3.1 Spot Amplitudes

In addition to computing the SDR for all the stars, we also examined
photometric variability measures. There have been previous suggestions of a
relationship between variability and rotation (e.g. Reinhold et al. (2013);
McQuillan, Mazeh & Aigrain (2014)); the results here are consistent with
previous analyses. Several variability diagnostics have been proposed and used
with Kepler data. They can be defined on various time scales, for example, or
with different methods of dealing with outliers. We employ the quarterly
“range” ($R_{var}$) that was used by Basri et al. (2011) and others as our
measure of variability. It is simply the difference in a quarterly light curve
between its 5% and 95% levels (after sorting all points by relative
intensity). Another popular choice is $S_{ph}$, which is the standard
deviation of the intensity over a few rotation periods (e.g. Mathur et al.
(2014)). They produce the same qualitative results; $R_{var}$ has an advantage
in general (though not here) that it does not depend on knowing the rotation
period. We computed $R_{var}$ (in ppt) for each quarter for a star, then used
the logarithm of the median quarterly $R_{var}$ value (for comparison with
SDR, which is also a logarithmic quantity) as our final measure of variability
for each star. The main difference with previous work is that we take a more
detailed look at what is happening for different stellar temperatures.

The $R_{var}$-Period relations are shown in Figs. 7 through 9. For the warmest
groups there is more scatter in the $R_{var}$-Period relationship, especially
at shorter periods (they also span the shortest set of rotation periods). This
is probably related to the beginnings of the tendency for stars with very
shallow convective envelopes to remain more rapidly rotating due to the lack
of magnetic braking, but they still exhibit some spot variability.

Figure 7: The relation between the rotation period (days) and the logarithm of
the median quarterly $R_{var}$ (in ppt) for the warmest temperature groups.
Note the changing period ranges in the different panels. An iterative cubic
fit is shown in red (with 2.5$\sigma$ outliers shown in purple). There is a
clear flattening at shorter periods for most temperature groups. A linear fit
to the longer half of the period range is also shown in green for each group.
The quiet and active Sun is indicated with the yellow diamonds.  Figure 8:
Same as Fig. 7 for the middle temperature groups.  Figure 9: Same as Fig. 7
for the cooler temperature groups. As with the SDR, there is a switch to
steeper upward slopes approaching short periods as one goes cooler.

As before we fit each group with a cubic polynomial in two iterations; in the
second we discarded the points that were more than 2.5$\sigma$ away from the
first fits. There is a tendency for the fits to flatten at rotation periods
shorter than 15-20 days. The flattening moves longward for the cooler stars,
which also extend to longer rotation periods. Between 4200-5800K the relation
for the more rapid rotators resembles a “saturation” behavior (very shallow
slope), as is seen for a number of more traditional activity diagnostics, such
as coronal X-rays (e.g. Reiners, Schüssler & Passegger (2014)). For cooler
stars (below 4200K) there is a group of more variable stars at periods less
than 10 days, which give the short period power law a steeper slope (the
opposite of what happens for the hotter stars), and the power law for these
stars at longer periods is shallower than for the hotter stars. These trends
continue down to 3200K. An odd dearth of lower values of $R_{var}$ appears
just before the half-way point in period between temperatures of 4200-5200K.
We are uncertain about the significance or meaning of that.

In order to compare with the linear fits for SDR-Period, we also produced
linear fits for the longer period half of each temperature group. Note that
the slopes of these longer period segments also show systematic behavior with
stellar temperature (quite similar to the SDR slopes in Fig. 6). This likely
means that $R_{var}$ and SDR are related, as we will argue more explicitly
below.

We also evaluated how the scatter in the $R_{var}$-Period relations compares
with the scatter in the SDR-Period relations. In Table 1 we tabulate the
standard deviations ($\sigma_{SDR}$,$\sigma_{Rng}$) for both SDR and $R_{var}$
for each temperature group after subtracting the cubic fits versus rotation
period. They are also shown as error bars in Fig. 6. The scatter for the stars
at the ends of the temperature ranges are highest, and the scatters for the
$R_{var}$ relations are slightly lower than for the SDR relations in the
middle temperature bins. Most of the relations have scatters of about a factor
of two (0.3 dex), but there is not an important difference in the dispersion
of the relations between SDR or $R_{var}$ versus period. This means that both
the SDR and $R_{var}$ are “rotation-activity” metrics, with similar levels of
efficacy. The SDR has more to do with the geometry of the starspot
distribution while $R_{var}$ has more to do with changes in starspot coverage,
although these are not completely independent of each other.

Table 1: Temperature Groups of Stars and Fit Parameters. Temp. Group (K) | Sample Size | Frac.(single)aaThe fraction of light curves in this temperature bin that are predominantly single mode.  | Slope (SDR) bbThe slopes of the linear fits to SDR in Figs.3 through 5.  | Slope ($R_{var}$) ccThe slopes of the longer period linear fits to $R_{var}$ in Figs.7 through 9.  |  $\sigma_{SDR}$ ddThe 1$\sigma$ residuals for the cubic fits to SDR in Figs. 3 through 5.  |  $\sigma_{Rng}$ eeThe 1$\sigma$ residuals for the cubic fits to $R_{var}$ in Figs.7 through 9.
---|---|---|---|---|---|---
3200-3400 | 158 | 0.316 | -0.027 | -0.015 | 0.374 | 0.412
3400-3600 | 456 | 0.281 | -0.026 | -0.014 | 0.373 | 0.312
3600-3800 | 673 | 0.279 | -0.031 | -0.010 | 0.356 | 0.286
3800-4000 | 877 | 0.239 | -0.033 | -0.015 | 0.335 | 0.291
4000-4200 | 1160 | 0.254 | -0.036 | -0.038 | 0.337 | 0.281
4200-4400 | 1652 | 0.245 | -0.035 | -0.039 | 0.333 | 0.274
4400-4600 | 983 | 0.255 | -0.039 | -0.041 | 0.316 | 0.282
4600-4800 | 1232 | 0.269 | -0.041 | -0.039 | 0.307 | 0.275
4800-5000 | 1807 | 0.316 | -0.045 | -0.047 | 0.308 | 0.264
5000-5200 | 3364 | 0.355 | -0.046 | -0.047 | 0.303 | 0.265
5200-5400 | 3650 | 0.358 | -0.049 | -0.049 | 0.287 | 0.271
5400-5600 | 3532 | 0.363 | -0.053 | -0.054 | 0.289 | 0.272
5600-5800 | 2877 | 0.342 | -0.054 | -0.055 | 0.298 | 0.310
5800-6000 | 2595 | 0.317 | -0.047 | -0.043 | 0.333 | 0.390
6000-6200 | 1613 | 0.274 | -0.029 | -0.046 | 0.360 | 0.407

Finally, we examine the amplitude of variability ($R_{var}$) in the single dip
segments compared with that in the double dip segments. Each ratio is computed
by taking the median $R_{var}$ of all the single segments and the median
$R_{var}$ of all the double segments in a given light curve. This is of
interest in part because two different scenarios for why double dip segments
appear could yield two different outcomes. If double dips arise primarily
because of changes in hemispheric spot coverage, then one might expect
$R_{var}$ to be somewhat similar in double dip and single dip segments. This
could happen either through additional spot coverage on the hemisphere that
was brighter during single dip times, or through reduction of the spot
coverage on the darker hemisphere to levels more like the (formerly) brighter
hemisphere. If, on the other hand, double dips arise because some of the spots
on the darker hemisphere during single dip times migrate over to the brighter
hemisphere (because of differential rotation), then the amplitudes of the
double dip segments should be systematically smaller than for single dip
segments. This would also yield the effect that the total spot coverage
(integrated light deficit) per rotation would stay relatively constant across
both single and double dip segments, whereas the first scenario would predict
changing coverage.

Figure 10: The ratios between the median $R_{var}$ of single segments to the
median $R_{var}$ of double segments for all the temperature groups. The curves
are smoothed versions of normalized histograms of the $R_{var}$ ratios for
each temperature group. The colors vary between reddest for the coolest group
to bluest for the warmest group. There is a strong tendency at most
temperatures for the ratio of single to double $R_{var}$ to be nearly twice as
high. More extreme ratios become increasingly common for cooler stars.

Figure 10 shows that there is a very clear preference for the single dip
median $R_{var}$ to be nearly twice as large as the double dip median
$R_{var}$. There are very few cases where they are equal or reversed. This
result is consistent with the fact that predominantly double dip stars also
have generally smaller median total $R_{var}$ (both SDR and median total
$R_{var}$ decrease with increasing period). The fact that changing the mix of
single/double segments will also change the mix of $R_{var}$ is probably the
explanation for the correlated behavior of the slopes in Fig. 6. There is a
significant fraction of single/double $R_{var}$ ratios that are even higher
than a factor of two for the coolest stars (higher contrast between single and
double dips), and the distribution of ratios is also a function of stellar
temperature as shown in Fig. 10. For stars above 5000K the single/double
$R_{var}$ ratios are more tightly clustered around a factor of slightly under
two. The tight clustering means that the ratio is maintained along the range
of values of total $R_{var}$; it does not matter much whether the total
variability is large or small. The broader distributions for the cooler groups
are driven by a tendency for the single $R_{var}$ to scatter to larger
relative values when the total variability is smaller, while the double
$R_{var}$ tends to continue to track well with the total $R_{var}$.

The ratio of a factor of two between the single and double $R_{var}$ is
roughly what would be expected if spots were simply redistributing themselves
on the star to generate either single or double dips, with relatively constant
coverage over a whole rotation (as argued above). Unfortunately the true
physical situation is not yet clear. As we discuss below, the preponderance of
higher amplitude single dips may instead (or in addition) be a result of a
tendency for larger numbers of spots to more likely result in a single dip
mode.

It is also not currently clear how the light deficit signal per rotation (true
total spot coverage) is actually changing, because there are still questions
about how the Kepler reductions deal with absolute photometric changes. The
pipeline differential light curves tend to have constant medians over each
quarter, and if one supposes that the unspotted intensity lies at some
constant value above the differential variations, that would imply relatively
unchanging spot coverage (and that the variations are largely due to spot
redistribution). But this may be an illusion due to the fact that real
absolute variations on timescales of a few weeks or more are suppressed by the
pipeline, in which case it is hard to know how the actual total light deficit
due to spots is changing (not to mention the possible influence of faculae on
this). We will study this question more in an upcoming paper that makes use of
absolute calibrations by Montet, Guadalupe & Foreman-Mackey (2017) for some
Kepler stars. For now all we can say for sure is that the differential
amplitudes of single dip segments tend to be about twice as large as those for
double dip segments.

###  3.2 Possible interpretations

One of the first questions these results raise is what sort of spot
configurations lead to the single dip mode or the double dip mode in our light
curves? For example, it is obvious that a configuration of 2 spots on opposite
sides of the star will usually lead to a double dip. Less obvious is what
minimum angular separation between the spots leads to a double dip. We
utilized the analytic spot modeling described by Walkowicz, Basri & Valenti
(2013) to begin to answer this. We started with $i=90^{\circ}$ models with 2
equatorial spots of equal size. In that case, the spots began to produce a
double dip (as detected by our procedure) when separated by more than
$105^{\circ}$. At a separation of $115^{\circ}$, however, the $i=30^{\circ}$
case is still classified as single, while higher inclinations are double. This
is an effect of the decreased contrast at lower inclination, and shows that
there are various complications in interpreting light curves. Changing the
ratio of the spot sizes or their latitudes provides further complications.

Before going further, it is appropriate to try to relate our results to those
of two recent papers which have some similarities in what they investigate.
The conclusions of Arkhypov et al. (2018) using about 1000 stars with the fast
half of our period range seem to boil down to an observation that the size of
the double dips (expressed through their half-period harmonic amplitude,
somewhat akin to the behavior of our double dip $R_{var}$ over time) are less
stable in time than the size of the single dips. We agree with that, but
believe it is premature to tie that directly to the behavior of sub-surface
diffusion and turbulence as they do. As we argue below, it is simply harder to
produce a double dip morphology in the light curve (requiring more special
spot configurations), and this naturally means that the double dips are more
variable in their amplitude. Furthermore, differential rotation by itself can
easily produce more rapidly changing double dips (in timing, depth, and
amplitude) and less rapidly changing single dips. That is because the double
dip amplitudes and positions are more sensitive to the actual spot
configuration while the single dips only require a predominance of one
hemisphere in coverage. As to the other conclusion of Arkhypov et al. (2018)
that cool slow rotators show more stable (in amplitude) single dips than
hotter faster rotators, we did not test this directly, but note that single
dips (which we agree are more stable) are definitely more common in faster
rotators at all temperatures.

Giles, Collier Cameron & Haywood (2017) worked with a similar number of stars
with similar temperatures, concentrated on periods of 10 and 20 days. They
also reference double dip segments, referring to them as “interpulse”
features. Because they perform auto-correlations of a whole light curve, the
signature of double dips is averaged in with single dips (and will be smaller
with larger SDR). Their primary interest is in how the auto-correlation
function decays as one looks a larger number of rotation periods out. Of
course, if the star is switching from single to double mode (or vice versa)
that will cause the auto-correlation function to decay faster. Given that the
double dip separations are more spread than the single dip separations (Fig.
2), stars with lower SDR will also degrade the auto-correlation function more
quickly. The primary conclusions of Giles, Collier Cameron & Haywood (2017)
are 1) the spot lifetime is longer for stars with greater variability, and 2)
cooler stars tend to have longer spot lifetimes (which is in good agreement
with Arkhypov et al. (2018)). In light of the our results, we could translate
those to 1) stars with larger $R_{var}$ spend more time in single mode, and 2)
cool stars tend to have larger $R_{var}$ and SDR at a given period. Both of
these conditions will produce less degradation of the auto-correlation
function. The first part of (2) was pointed out by Basri et al. (2011) and
reinforced by several later papers. The second part is partially true; the
single fraction does increase for the coolest groups compared to stars of
middle temperatures (Fig. 6), although stars near solar temperature have even
larger single fractions (although these are over the full period range).
Perhaps it is the combination of the single fraction and $R_{var}$ that
produces the cited result for the decay of the auto-correlation function for
cool stars.

We do not mean to imply that the conclusions of Giles, Collier Cameron &
Haywood (2017) are wrong (indeed, we agree below that spot lifetimes may well
be longer for higher SDR stars). We just want to point out that it is not yet
clear whether spot lifetimes or spot configurations are being tested. Many
authors have interpreted light curves in simple terms, imagining that the
presence or absence of a secondary dip at about half a period reflects the
behavior of a secondary spot on the other side of the star from a primary
spot. We argue that this is overly simplistic, and can easily miss other
different physical effects such as migration of some spots from one side to
the other. In that case it would be a mistake to infer spot lifetimes from the
behavior of the light curve. It is of interest to study in more detail how
long each single/double mode lasts, and that would be one good thing to study
next. Clearly stars that have an SDR of 0.5 or more spend most of their time
in single mode, so their typical single segment must be longer than for lower
SDR values.

To properly understand starspot light curve morphologies, in a following paper
we will describe a far more complete exploration of the many-parameter space
that spot models can occupy. Among the most relevant parameters are 1) basic
spot characteristics like number, size distribution, and contrast, 2) the
distributions in longitude and latitude: over what range?, are they random or
grouped?, are they confined to stripes or belts?, 3) spot evolution: how often
do spots appear? what are the growth and decay timescales?, 4) differential
rotation: do they have different rotation periods? by how much?, how does it
depend on their latitude?, and 5) the stellar inclination. Considerations such
as these have not been carefully considered by most authors; one would like to
know how they are manifested in metrics in intensity and frequency space.
Metrics of the light curve morphology that should be studied include some
measure of differential amplitude (like $R_{var}$), inferred spot coverage
(which requires also knowing the unspotted “continuum”), the SDR, periodogram
structures (like double peaks), other measures of complexity (like auto-
correlation functions and dip separation and depth histograms) and
evolutionary behavior of the metrics. For each parameter set, we run thousands
of trials with certain parameters randomized, and look at the statistical
behavior of the metrics. Those results are far too lengthy to include in this
paper, but we next present some preliminary information on how the SDR
behaves.

One of the main conclusions from our extensive spot model trials is that the
more common result in a general exploration of parameter space is a preference
for positive SDR values. In other words, single dip light curve segments are
decidedly more common than double dip segments for most of the configurations
we have tried. This is in tension with the empirical results presented here
(third column of Table 1), which show that most stars (other than the shortest
period ones) tend to favor double dip morphologies. The only part of parameter
space we have found so far that clearly favors double dips requires spot
evolution with only a few spot groups present at a time and spot lifetimes not
much over 1-2 rotation periods. As it happens, this describes the current
state of the Sun (which indeed presents a strongly double dip light curve
structure; Fig. 3). A sensible inference is that the strongly double dip stars
(generally the longer period cases) have few spot groups, which probably don’t
last too long. Conversely, the fact that short period stars tend to have
single dip structures and larger $R_{var}$ is consistent with the proposition
that they have larger spot coverage (many spots present at the same time)
and/or longer spot lifetimes. A remaining unsettled question is whether the
spot coverage is relatively constant or changes substantially as a given star
switches back and forth between single to double mode (resolving this requires
a better understanding of how to interpret Kepler light curves, which we are
also working on).

Our tentative overall interpretation of our results is that the relation
between SDR and rotation may be a diagnostic of a combination of spot coverage
and lifetime, with both decreasing as the rotation period increases. Such a
correlation is not surprising in a qualitative sense, but the combination of
coverage and lifetime somehow manages to be closely related to the stellar
rotation period. We do not yet understand why there are clear relationships
between stellar surface temperature (mass) and the steepness of the SDR-Period
relations. It is tempting to suppose that the Rossby number is somehow
involved, since that is a widely used variable that introduces the stellar
temperature (through the convective overturn time) into considerations of
magnetic activity. We also are not currently certain we know why the single
segment $R_{var}$ is commonly nearly twice the double segment $R_{var}$ within
a light curve. These mysteries suggest that one might learn something about
the stellar magnetic dynamo, or at least its production of starspots, from a
better understanding of these relationships. Starspots may probe the
connections between the depth (and perhaps vigor) of the convection zone and
the geometry of field production. Spot lifetimes do indeed depend on diffusive
processes that are related to sub-surface motions. That is the suggestion in
the work of Arkhypov et al. (2018), and although we are not currently
convinced that their explanations have already been demonstrated, they are
heading in the right direction. Clearly we need to better understand what the
single/double nature of starspot light curves is really telling us. This paper
is an early step in the analysis of that subject.

This paper includes data collected by the Kepler mission. Funding for the
Kepler mission was provided by the NASA Science Mission directorate. Most of
the data presented in this paper were obtained from the Mikulski Archive for
Space Telescopes (MAST). STScI is operated by the Association of Universities
for Research in Astronomy, Inc., under NASA contract NAS5-26555. This research
has also made use of the NASA Exoplanet Archive, which is operated by the
California Institute of Technology, under contract with the National
Aeronautics and Space Administration under the Exoplanet Exploration Program.
We would like to thank the referee for a careful reading of the manuscript and
a number of suggestions that improved the paper. Facilities: Kepler, MAST,
Exoplanet Archive.

## References

  * Akeson et al. (2013) Akeson, R.L., Chen, X., Ciardi, D., Crane, M., Good, J., 2013, PASP, 125, 898
  * Arkhypov et al. (2018) Arkhypov, O.V., Khodachenko, M.L., Lammer, H., Güdel, M., Lüftinger, T., Johnstone, C.P., 2018, MNRAS, 473, L84
  * Baglin (2003) Baglin, A., 2003, Astrophys. Space Phys. Res., 31, 345
  * Basri et al. (2011) Basri, G., Walkowicz, L., Batalha, N., Gilliland, R. Jenkins, J., Borucki, W., Koch, D., et al. 2011, ApJ, 713, L155
  * Borucki et al. (2010) Borucki, W.J., Koch, D., Basri, G., Batalha, N., Brown, T. et al., 2010, Science, 327, 977
  * Giles, Collier Cameron & Haywood (2017) Giles, H.A.C., Collier Cameron, A., Haywood, R.D., 2017, MNRAS, 472, 1618
  * Mathur et al. (2014) Mathur, S., García, R.A., Ballot, J., Ceillier, T., Salabert, D., 2014, A&A, 562, 124
  * McQuillan, Mazeh & Aigrain (2014) McQuillan, A., Mazeh, T., Aigrain, S., 2014, ApJS, 211, 24
  * Meibom et al. (2015) Meibom, S., Barnes, S.A., Platais, I., Gilliland, R.L., Latham, D.W., Mathieu, R.D., 2015, Nature, 517, 589
  * Montet, Guadalupe & Foreman-Mackey (2017) Montet, B. T., Guadalupe, T., Foreman-Mackey, D., 2017,ApJ, 851, 116
  * Noyes et al. (1984) Noyes, R.W., Hartmann, L.W., Baliunas, S.L., Duncan, D.K., Vaughan, A.H., 1984. ApJ, 279, 763
  * Reiners, Schüssler & Passegger (2014) Reiners, A., Schüssler, M., Passegger, V.M., 2014, ApJ, 794, 144
  * Reinhold et al. (2013) Reinhold, T., Reiners, A., Basri, G., 2013, A&A, 560, 4
  * See et al. (2017) See, V., Jardine, M., Vidotto, A.A., Donati, J.-F., Boro Saikia, S., Fares, R., 2017, MNRAS, 466, 1542
  * Semel (1989) Semel, M., 1989, A&A, 225, 456
  * Strassmeier (2009) Strassmeier, K., 2009, ARA&A, 17, 251
  * Vogt & Penrod (1983) Vogt, S.S., Penrod, G.D., 1983, PASP, 94, 565
  * Walkowicz, Basri & Valenti (2013) Walkowicz, L., Basri, G., Valenti, J., 2013, ApJS, 205, 17


11institutetext: Leibniz-Institut für Astrophysik Potsdam (AIP), An der
Sternwarte 16, 14482 Potsdam, Germany
11email<EMAIL_ADDRESS>22institutetext: Institut für Physik und Astronomie,
Universität Potsdam, Karl-Liebknecht-Str. 24/25, 14476 Potsdam, Germany
33institutetext: Zentrum für Astronomie der Universität Heidelberg,
Landessternwarte, Königstuhl 12, 69117 Heidelberg, Germany 44institutetext:
Max Planck Institute for Astronomy, Königstuhl 17, 69117, Heidelberg, Germany
55institutetext: Instituto de Astronomía, Universidad Nacional Autónoma de
México, A. P. 106, C.P. 22800 Ensenada, B. C., México 66institutetext:
Dipartimento di Fisica e Astronomia, Università degli Studi di Bologna, Via
Gobetti 93/2, I-40129 Bologna, Italy 77institutetext: Freie Universität
Berlin, Fachbereich Physik, Arnimallee 14, 14195 Berlin, Germany

# Insights from Super-Metal-Rich Stars:
Is the Milky Way bar young?

S. Nepal 11 2 2 C. Chiappini 11 G. Guiglion 3344 M. Steinmetz 11 A. Pérez-
Villegas 55
A. B. Queiroz 11 A. Miglio 66 P. Dohme 77 A. Khalatyan 11

(Received 23 October 2023 / accepted 28 November 2023)

###### Abstract

Context. Bar formation and merger events can contribute to the rearrangement
of stars within the Galaxy in addition to triggering star formation (SF)
epochs. Super-metal-rich (SMR) stars found at the Solar Neighbourhood (SNd)
can be used as tracers of such events as they are expected to originate only
in the inner Galaxy and have definitely migrated.

Aims. We aim at studying a homogeneous and large sample of SMR stars in the
SNd to provide tighter constraints on the epoch of the bar formation and its
impact on the Milky Way (MW) disc stellar populations.

Methods. We investigate a sample of 169 701 main sequence turnoff (MSTO) and
sub-giant branch (SGB) stars with 6D phase space information and high-quality
stellar parameters coming from the hybrid-CNN analysis of the _Gaia_ -DR3 RVS
stars. We compute distances and ages using the StarHorse code with a mean
precision of 1% and 11%, respectively. From these, 11 848 stars have
metallicity ([Fe/H]) above 0.15 dex.

Results. We report a metallicity dependence of spatial distribution of stellar
orbits shown by the bimodal distribution in the guiding radius
($\mathrm{R_{g}}$) at 6.9 and 7.9 kpc, first appearing at [Fe/H]$\sim$0.1 dex,
becoming very pronounced at larger [Fe/H]. In addition, we’ve observed a trend
where the most metal-rich stars, with [Fe/H]$\sim$0.4 dex, are predominantly
old (9-12 Gyrs) but show a gradual decline in [Fe/H] with age, reaching around
0.25 dex at about 4 Gyrs ago, followed by a sharp drop around 3 Gyrs ago.
Furthermore, our full dataset reveals a clear peak in the age-metallicity
relationship during the same period, indicating a SF burst around 3-4 Gyrs ago
with slightly sub-solar [Fe/H] and enhanced [$\alpha$/Fe].

Conclusions. We show the SMR stars are good tracers of the bar activity. We
interpret the steep decrease in number of SMR stars at around 3 Gyr as the end
of the bar formation epoch. In this scenario, the peak of bar activity also
coincides with a peak in the SF activity in the disc. Although the SF burst
around 3 Gyr ago has been reported previously, its origin was unclear. Here,
we suggest the SF burst to have been triggered by the high bar activity, 3-4
Gyr ago. According to these results and interpretation, the MW bar could be
young.

###### Key Words.:

Galaxy: abundances - Galaxy: evolution - Galaxy: kinematics and dynamics -
stars: fundamental parameters

##  1 Introduction

Stars are the luminous storytellers of our Galactic saga. By studying the
properties of these stars, their ages, chemical abundances and motions, we can
trace back the history of our Galaxy (Pagel 1997; Matteucci 2001, 2021;
Freeman & Bland-Hawthorn 2002). Among the different stellar populations in the
Milky Way (MW), super-metal-rich (SMR111SMR stars are defined such that their
metal-abundance exceeds the metallicity of the local present-day interstellar
medium, and because of the negative radial abundance gradient in the MW, its
value varies with galactocentric distance.) stars are very interesting as they
are expected to be formed only in the inner regions of our Galaxy out of
materials enriched by previous generations of stars (Grenon 1972; Steinmetz &
Mueller 1994; Trevisan et al. 2011; Casagrande et al. 2011; Kordopatis et al.
2015; Anders et al. 2017; Miglio et al. 2021). These SMR stars currently
residing in the SNd have definitely migrated from the inner Galaxy (Feuillet
et al. 2018; Chen et al. 2019; Dantas et al. 2023) and their study can inform
us about the processes that bring them to the outer disc. Indeed, Queiroz et
al. (2021) have shown the inner kpcs of the Galaxy to provide a large
reservoir of high-metallicity stars (either as part of an inner-thin disc or
on bar shape orbits), and these are primary good candidates for migration.

The Galactic bar is considered an important perturber for generating radial
motion of stars and gas in the galactic disc (e.g. Sellwood 2014). According
to N-body simulations, during the formation and the phase of strong bar
activity, stars in the disc are significantly redistributed, with the highest
probability of migration at the bar resonances (Minchev & Famaey 2010; Halle
et al. 2015; Khoperskov et al. 2020a). However, an open question still remains
on the epoch of MW bar formation.

Old to young formation times for the MW bar have been proposed in both
observational and simulation studies conducted so far. Bovy et al. (2019) and
Wylie et al. (2022), using red giant stars in the inner Galaxy with ages
estimated using astroNN (Leung & Bovy 2019), proposed an older bar at least 7
to 8 Gyr old based on the mean and peak of the age distribution. Cole &
Weinberg (2002), using Carbon stars that trace the bar, suggested that bar
likely formed $\sim$3 Gyr ago. However, it is important to recognise that the
bar is a cumulative entity, it will indiscriminately contain older stars and
stars born during bar formation. So, care has to be taken to differentiate
older stars from those formed during bar formation (e.g. see de Sá-Freitas et
al. 2023a). Tepper-Garcia et al. (2021), using a tailored N-body model of
barred MW, suggested that the bar formed 3-4 Gyrs ago. In context of external
galaxies, de Sá-Freitas et al. (2023b) recently reported the discovery of
young bars, formed 4.5 and 0.7 Gyr ago, for spiral galaxies NGC 289 and NGC
1566 using the SF history of nuclear discs. This suggests that some disc
galaxies, with stellar mass comparable to MW, could settle on longer time-
scales.

However, because SMR stars are rare, we still lack a statistically significant
sample to constrain this process. This situation has radically changed thanks
to third data release (DR3) of the ESA _Gaia_ mission (Gaia Collaboration et
al. 2023). _Gaia_ -DR3 has provided about one million spectra from the Radial
Velocity Spectrometer (RVS)222The RVS spectra were originally analysed during
_Gaia_ DR3 (10.17876/gaia/dr.3) by the General Stellar Parametriser for
spectroscopy (GSP-Spec, Recio-Blanco et al. 2023) module of the Astrophysical
parameters inference system (Apsis, Creevey et al. 2023). of which only
$\sim$178,000 stars have the good quality stellar parameters after the
applying the recommended flags (Recio-Blanco et al. 2023). A large portion of
the published RVS spectra are of low signal-to-noise (S/N; 15 - 25 per pixel)
and can be challenging for traditional spectroscopic pipelines. Guiglion et
al. (2023, G23) using a hybrid Convolutional Neural-Network (CNN) method,
reanalyzed the RVS sample to derive atmospheric parameters
($\mathrm{\textit{T}_{\text{eff}}}$, log(g), and overall [M/H]) and chemical
abundances ([Fe/H] and [$\alpha$/M]), by supplementing extra information from
_Gaia_ magnitudes (Riello et al. 2021), parallaxes (Lindegren et al. 2021) and
XP coefficients (De Angeli et al. 2023). The hybrid-CNN was trained with a
high-quality training sample based on APOGEE DR17 (Abdurro’uf et al. 2022)
labels and show precision and accuracy comparable to external data-sets such
as GALAH and asteroseismology. Thanks to the novel method, G23 significantly
improved the number of reliable targets that can be used for Galactic
archaeology. The G23 catalog, has a large number of SMR stars, including a
sample of MSTO and SGB stars, enabling the use of these traces to constrain
the epoch of MW bar formation.

In this Letter, we explore the age-[Fe/H] (AMR) and the age-[$\alpha$/Fe]
(AAR) relationships along with the chemo-dynamics for a sample of MSTO and SGB
stars selected from the RVS-CNN catalog of G23. In Sect. 2, we describe our
sample focusing on the methods used to obtain the stellar ages and kinematics.
In Sect. 3 we present our results and in Sect. 4 we present the main
conclusions.

##  2 Data

We obtain a sample of stars, with $-$0.7 ¡ [Fe/H]¡ 0.5 including a large set
of the SMR stars, from the catalog of G23. We selected MSTO and SGB stars (see
Queiroz et al. 2023) after the application of the G23 recommended flags. We
keep only the MSTO-SGB stars as the stellar ages from isochrone fitting
methods are most reliable for these evolutionary stages (e.g. Soderblom 2010).
We computed the distances and stellar ages with the StarHorse Bayesian
isochrone-fitting method (Queiroz et al. 2018, 2023) and integrate the orbits
of the stars using Galpy (Bovy 2015). For details on these computations see
Appendix A.

We selected stars with relative age uncertainty less than 20%, distance
uncertainty less than 5% and extinction uncertainty less than 0.2 Mag. We are
left with a sample with mean uncertainty of 11% for age and 1% for distance
due to very low parallax errors in the extended SNd thanks to _Gaia_. We
further remove any stars with poor astrometric solutions by limiting RUWE ¡
1.4 and also remove known variable stars by using _Gaia_ flag
’phot_variable_flag’$\neq$’VARIABLE’ (see Gaia Collaboration et al. 2023).

Figure 1: The properties of our selected sample - a) Distance from the
galactic mid-plane (Z) vs the Galactocentric distance (R); b) the Kiel diagram
(log(g) vs $\mathrm{\textit{T}_{\text{eff}}}$); and c) the [$\alpha$/Fe] vs
[Fe/H] diagram for the sample, colored by number density.

This gives us a final sample of 169 727 stars, including 19 367 stars with
[Fe/H]¿0.1 dex. As a test, we apply the 13 recommended flags on the GSP-Spec
parameters (Recio-Blanco et al. 2023) to find a total sample reduced to 20 269
stars with only 4 853 stars with [Fe/H]¿0.1 dex.

In Fig.1, we present our sample properties. Although the stars are widely
distributed in Z vs R space, the number of stars decreases as we move away
from SNd, as expected for a sample of MSTO+SGB stars (see Queiroz et al.
2023). Therefore the sample is essentially dominated by thin disc stars
(low-[$\alpha$/Fe] population).

Figure 2: Distribution of the guiding radii ($\mathrm{R_{g}}$) in bins of
metallicity ([Fe/H]). The [Fe/H] values increase from top left to bottom
right, i.e. from panel a) to panel i) (number of stars and the [Fe/H] range
shown). The two dotted lines at 6.9 and 7.9 kpc represent the two peaks of the
$\mathrm{R_{g}}$ distribution.

##  3 Results and Discussion

###  3.1 Super-metal-rich stars as tracers of MW bar activity

Figure 2 shows the distribution of guiding radius ($\mathrm{R_{g}}$), in bins
of metallicity, for stars in the range $-$0.3¡[Fe/H]¡0.5. The [Fe/H] values
increase as we move from panel a) to panel i). These histograms show that a
clear bimodality in the distribution of $\mathrm{R_{g}}$ appears as we move
from sub-solar [Fe/H] to SMR stars. The bimodality visibly appear first in the
bin 0.10¡[Fe/H]¡0.15, as seen in panel (d), and becomes stronger for the
higher [Fe/H] bins. An inner and an outer peak is seen at 6.9 and 7.9 kpc, for
all super-solar bins. We validate this finding with an external sample of red-
giant stars in Appendix D. We note that the two peaks are not equally
populated, with higher number of stars near the solar radius, and attribute
this to the _Gaia_ selection function.

Unlike SMR stars, sub-solar to solar [Fe/H] stars are equally plausible to be
formed in the SNd as well as the inner and outer galaxy, so their actual birth
location can be uncertain due to radial migration (e.g. Minchev et al. 2013,
2014 and references therein). We also know the innermost regions of the MW to
be populated of very metal-rich stars (e.g see Barbuy et al. 2018 and
references therein). In particular, Queiroz et al. (2021), using APOGEE and
_Gaia_ surveys to study bar and bulge region, found a large cache of metal-
rich stars in inner Galaxy. These metal-rich stars, found in cold to highly
eccentric orbits, could be a source of origin for our local SMR stars.

What could be the mechanism that brought them here? The SMR stars in our
sample, most of them in cold orbits with 96% with ecc¡0.3, distributed
favourably at certain $\mathrm{R_{g}}$ entices us to view their distribution
in context of evolution and dynamical process in the Galaxy (See Appendix B
for a comparison of velocity dispersion of metal-rich stars with the full
sample.). The Galactic bar has been long considered an efficient source of
radial migration (Combes & Sanders 1981; Pfenniger & Norman 1990). Galaxy
models and N-body simulations have shown that slow down of Galactic bars can
lead to migration of stars from inner to outer galaxy by trapping them in
resonances as they travel outwards in the disc (e.g. Athanassoula 2003;
Khoperskov et al. 2020a; Chiba & Schönrich 2021).

Although dependent on the exact dynamical recipe, many recent studies have
placed bar Corotation (CR) between radius of 6-7 kpc (e.g. Portail et al.
2017; Khoperskov et al. 2020b; Chiba & Schönrich 2021) and the Local arm at
around 8 kpc (e.g. Palicio et al. 2023) linking it to Outer Lindbald
Resonances (OLR; Fragkoudi et al. 2019; Khoperskov et al. 2020b). Furthermore,
Chen et al. (2022) using SMR stars, find ridges and undulations in the $\phi$
vs $L_{Z}$ plane similar to the orbits trapped in resonances of a slow bar as
in the model of Monari et al. (2019). See Appendix F for test on how metal-
rich stars are better suitable for probing bar activity. Therefore, the inner
and outer $\mathrm{R_{g}}$ peaks, observed for the SMR stars at 6.9 and 7.9
kpc, most probably correspond to the imprints of the Galactic bar.

From our observations (see Fig.2) we conclude that the SMR stars, currently in
the SNd, trace the signatures of bar activity. In the next subsection we
employ these SMR stars, with good ages, to study the formation epoch of the MW
bar.

Figure 3: Left panel: Stellar ages as a function of guiding radii for stars
with [Fe/H]¿0.1 dex. The colors represent number of stars per bin. A Kernel
Density Estimate (KDE) has also been over-plotted. A dashed black line at 2.7
Gyr is also shown. Right panel: Upper envelope of the age-metallicity relation
traced by our full sample. The triangles represent the median [Fe/H] of 20
most metal-rich stars in each age bin and the error-bars show the scatter in
[Fe/H].

###  3.2 Timing the bar formation epoch with the youngest super-metal-rich
stars

In the left panel of Fig. 3 we present the stellar ages as a function of
$\mathrm{R_{g}}$, for the 19 367 metal-rich stars for which bimodality is
clearly seen in Fig. 2, i.e. for [Fe/H]¿0.1 dex. The plot shows that the
metal-rich stars have a wide range of ages with higher prevalence at older
ages (11 to 6 Gyr). This hints to the fact that a large number of metal-rich
stars were formed at early times in our Galaxy’s history. The bimodal nature
of $\mathrm{R_{g}}$, with peaks at 6.9 and 7.9 kpc, is clearly seen for all
ages and the spread in $\mathrm{R_{g}}$ is narrower for the younger stars.
Although some of these older stars could have already dispersed to the outer
galaxy at earlier times through various processes such as merger events (e.g.
Helmi 2020), most of them are still in cold orbits. Interestingly, we see a
near absence of metal-rich stars with ages younger than $\sim$2.7 Gyr old (see
the dashed line in Fig. 3, also see Fig. 9). This cannot be attributed to
selection effects on our sample such as color or magnitude limits, see
Appendix C. Also, significant number younger stars with sub-solar [Fe/H] are
present in our sample, see Figs. 4 and 9. Therefore, the dearth of younger SMR
stars could be related to the mechanism that brought these stars formed in the
inner Galaxy to the SNd.

In the right panel of Fig. 3 we present the time evolution of the upper
envelope of the metallicity for our whole sample. The plot reveals that our
Galaxy quickly attained the highest levels of metallicity enrichment,
levelling around [Fe/H]=0.4 dex, already at $\sim$10 Gyr ago. Then starting at
8.5 Gyr ago we observe a slow decreases, at $-$0.04 dex/Gyr, up to 4 Gyr. Then
for a period of 1 Gyr, from 4 to 3 Gyr ago, we see a plateau at [Fe/H]=0.25
dex, highlighted by an orange band. After 3 Gyrs we see a steep decline in the
[Fe/H] envelope at a rate of $-$0.15 dex/Gyr (see Appendix D for a validation
with an asteroseismic sample). This decline of the [Fe/H] envelope further
confirms that the mechanism, that brought majority of these metal-rich stars
to SNd, has ceased. We note here that Minchev et al. (2011) identified a
radial migration mechanism, via bar and spiral structure interaction, that can
migrate stars throughout the disc of MW-like galaxy in ¡ 1 Gyr period.
However, the migration induced by the bar is most effective at the epoch of
bar growth and, therefore, stars in the inner Galaxy can be moved en masse to
larger $\mathrm{R_{g}}$ (e.g. Halle et al. 2015).

Considering MW bar growth and outward propagation of the resonances as the
transport mechanism of these SMR stars to the SNd (Minchev & Famaey 2010;
Khoperskov et al. 2020a; Iles et al. 2023), we can use the stellar ages of
these youngest SMR stars to place a constraint on the epoch of intense bar
growth following bar formation. It is worth noting that the MW bar, and inner
Galaxy in general, when currently observed will contain stars of all ages from
very old in higher number to very young, including those formed before, during
and after the epoch of bar formation. However, only the youngest metal-rich
stars, formed during the epoch of bar formation, can be confidently employed
to constrain the age of bar. Based on this absence of younger SMR stars, we
place the intense bar growth for the MW ending at 2.5 to 3 Gyrs ago.

Figure 4: Top (a): 2D density distribution showing [Fe/H] as a function of
stellar ages for the full sample. Bottom (b): 2D density distribution showing
[$\alpha$/Fe] as a function of stellar ages. A Kernel Density Estimate (KDE)
have also been overplotted to highlight the density features in the AMR and
AAR relationships. The colors represent number of stars per bin in log scale.

###  3.3 Enhanced star formation triggered by the bar formation

Considering that intense bar growth for MW ended at $\sim$3 Gyrs ago (see
previous Section), we can expect an epoch of bar enhanced SF succeeding this
event (e.g. Baba & Kawata 2020). In this section, we explore SF in the Age-
Metallicity and Age-Alpha Relationships (AMR and AAR) and look for possible
imprints of such an event.

Fig. 4 shows, the AMR (panel a) and AAR (panel b) stellar density
distributions for the whole sample. The AMR shows three distinct density
features, at correspondingly three age regimes.

At the oldest regime, from 12 to $\sim$9-10 Gyrs, an oblique feature shows a
steep rise in [Fe/H] with age with the broadest [Fe/H] distribution. This
oblique feature depicts an epoch of intense SF that leads to a rapid chemical
enrichment in the early MW from the lowest to the highest [Fe/H] , i.e. from
$-$0.7 to 0.4 dex. In the accompanying AAR, at this regime, we see the high
values of [$\alpha$/Fe] corresponding to the SF burst in the early galaxy
leading to the formation of the chemical thick disc and old bulge (e.g.
Fuhrmann 1998; Anders et al. 2018; Miglio et al. 2021; Montalbán et al. 2021;
Queiroz et al. 2023; Xiang & Rix 2022). In addition, this epoch is also
attributed to the merger of MW with Gaia-Sausage-Enceladus (GSE; Belokurov et
al. 2018; Helmi et al. 2018).

Between 9 to 5 Gyrs we see a horizontal blob with most of the stars between
$-$0.2¡[Fe/H]¡0.0. This epoch reflects a slow and steady SF and is attributed
to the growth of the thin disc (e.g. Chiappini et al. 1997; Minchev et al.
2013). This regime shows a birth of a large number of stars with sub-solar
[Fe/H]. A trend of decreasing [$\alpha$/Fe] with age is seen, also a sign of
slow but continued SF. Growth of the MW thin disc in an inside-out manner, due
to steady gas infall is expected in this regime (see Matteucci 2021 and
references therein).

Between 5 to 4 Gyrs we see a decline in SF, followed by another feature
between 4 to 2.7 Gyr. This oblique feature shows an enhanced SF with a slight
decrease in metallicity at younger ages. An increase in the [$\alpha$/Fe] is
also seen, showing a signature of SF burst. This result further confirms the
enhanced SF during similar epoch previously reported (Rocha-Pinto et al. 2000;
Isern 2019; Mor et al. 2019; Sysoliatina & Just 2021; Sahlholdt et al. 2022;
Imig et al. 2023). These mentioned previous works have not conclusively
attributed any physical process to this SF phase. Considering that intense bar
growth period ended at $\sim$3 Gyr ago, see Sect. 3.2, we deduce that this
star-burst is caused by the high activity of the bar evolution.

A likely scenario one could consider is that this high bar activity is simply
a buckling of an old bar. But, there is still no consensus if the bar-buckling
triggers SF in the disc (See Fragkoudi et al. 2020 and for opposite view see
Debattista et al. 2006; Łokas 2020).

Interestingly, in their study of the outer bar region, Wylie et al. (2022)
also find an abrupt decline of younger stars followed by a significant
fraction of stars at 2-4 Gyr for the inner disc and ring region (see their
Fig. 3). Additionally, bar formation triggered star burst have been seen to
typically last for a duration of $\sim$1 Gyr (e.g. see Baba & Kawata 2020),
this corroborates our assessment. Furthermore, bars have been observed to
enhance SF in external galaxies. For example, see study by Lin et al. (2020)
on low-redshift galaxies using integral field spectroscopy (also see Ellison
et al. 2011). Hence, we conclude that MW bar formation occurs at $\sim$4 Gyr
ago with end of strong bar activity at $\sim$2.7 Gyr ago.

##  4 Conclusion

We have explored a large and homogeneous sample of 169 701 MSTO and SGB stars
with 6D phase space information and high-quality stellar parameters coming
from the _Gaia_ -DR3 RVS analysis of Guiglion et al. (2023). We supplemented
the chemical abundances with stellar ages, distances and kinematics to study
the epoch of MW bar formation. Thanks to _Gaia_ DR3 we obtain a mean distance
uncertainty of 1% which greatly contributes to the findings of this Letter.

The new data shows two new results: a) a clear bimodality in $\mathrm{R_{g}}$
at all ages (¿3 Gyrs) of SMR stars, and b) a dearth of SMR stars younger than
3 Gyr. These results imply:

  * •

Milky Way’s bar had a strong activity phase lasting $\sim$1 Gyr ending at
$\sim$3 Gyr ago. During the phase of strong bar activity, stars formed in the
inner region (bulge/bar), are significantly redistributed across the outer
disc, with the highest probability of migration around bar resonances. We
verify this with a dearth of SMR stars younger than 3 Gyrs and an observed
bimodality in the $\mathrm{R_{g}}$ around the bar resonances for the SMR
stars, which are, with high confidence, formed in the inner galaxy and brought
here during the strong activity phase of bar formation.

  * •

We detect an enhancement in the global SF (around $-$0.3 ¡ [Fe/H]¡ 0.0) at
around 3 Gyr which, by that time had already declined in the local thin disc.
Although this SF enhancement has been detected previously in the literature,
the age coincidence with our estimate of the bar age suggest these events to
be related.

  * •

We suggest that, due to mixing and strong gas inflow due to bar there was an
epoch of enhanced SF which is seen in Fig. 4. This gas inflow during this
phase lowers the upper floor of [Fe/H] vs age distribution while causing an
increase in the [$\alpha$/Fe] due to intense SF.

Future spectroscopic surveys such as 4MIDABLE-LR Chiappini et al. (2019) will
enable a large increase in the number of super metal rich stars with full 6D
phase space and chemical information. In addition, the Japan Astrometry
Satellite Mission for INfrared Exploration (JASMINE - Kawata et al. 2023) will
be able to provide additional constraints on the age of the MW bar by
precisely measuring the age of the Nuclear Stellar Disc.

###### Acknowledgements.

SN thanks the E-science & IT team for COLAB service, computational clusters
and research infrastructure at AIP. SN thanks Sergey Khoperskov, Friedrich
Anders and Ivan Minchev for the helpful suggestions. GG acknowledges support
by Deutsche Forschungs-gemeinschaft (DFG, German Research Foundation) –
project-IDs: eBer-22-59652 (GU 2240/1-1). This project has also received
additional funding from the European Research Council (ERC) under the European
Union’s Horizon 2020 research and innovation programme (Grant agreement No.
949173). APV acknowledges the DGAPA–PAPIIT grant IA103122. AM acknowledges
support from the ERC Consolidator Grant funding scheme (project
ASTEROCHRONOMETRY, G.A. n. 772293). 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 work made use of
overleaf (https://www.overleaf.com/), and of the following python packages:
matplotlib (Hunter 2007), numpy (Harris et al. 2020), pandas (McKinney 2010),
seaborn (Waskom 2021). This work also benefited from topcat (Taylor 2005).

## References

  * Abdurro’uf et al. (2022) Abdurro’uf, Accetta, K., Aerts, C., et al. 2022, ApJS, 259, 35
  * Anders et al. (2017) Anders, F., Chiappini, C., Minchev, I., et al. 2017, A&A, 600, A70
  * Anders et al. (2018) Anders, F., Chiappini, C., Santiago, B. X., et al. 2018, A&A, 619, A125
  * Astropy Collaboration et al. (2022) Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, ApJ, 935, 167
  * Athanassoula (2003) Athanassoula, E. 2003, MNRAS, 341, 1179
  * Baba & Kawata (2020) Baba, J. & Kawata, D. 2020, MNRAS, 492, 4500
  * Barbuy et al. (2018) Barbuy, B., Chiappini, C., & Gerhard, O. 2018, ARA&A, 56, 223
  * Belokurov et al. (2018) Belokurov, V., Erkal, D., Evans, N. W., Koposov, S. E., & Deason, A. J. 2018, MNRAS, 478, 611
  * Bovy (2015) Bovy, J. 2015, ApJS, 216, 29
  * Bovy et al. (2019) Bovy, J., Leung, H. W., Hunt, J. A. S., et al. 2019, MNRAS, 490, 4740
  * Bressan et al. (2012) Bressan, A., Marigo, P., Girardi, L., et al. 2012, MNRAS, 427, 127
  * Buder et al. (2021) Buder, S., Sharma, S., Kos, J., et al. 2021, MNRAS, 506, 150
  * Casagrande et al. (2011) Casagrande, L., Schönrich, R., Asplund, M., et al. 2011, A&A, 530, A138
  * Chen et al. (2022) Chen, Y., Zhao, G., & Zhang, H. 2022, ApJ, 936, L7
  * Chen et al. (2019) Chen, Y. Q., Zhao, G., Zhao, J. K., et al. 2019, AJ, 158, 249
  * Chiappini et al. (1997) Chiappini, C., Matteucci, F., & Gratton, R. 1997, ApJ, 477, 765
  * Chiappini et al. (2019) Chiappini, C., Minchev, I., Starkenburg, E., et al. 2019, The Messenger, 175, 30
  * Chiba & Schönrich (2021) Chiba, R. & Schönrich, R. 2021, MNRAS, 505, 2412
  * Cole & Weinberg (2002) Cole, A. A. & Weinberg, M. D. 2002, ApJ, 574, L43
  * Combes & Sanders (1981) Combes, F. & Sanders, R. H. 1981, A&A, 96, 164
  * Creevey et al. (2023) Creevey, O. L., Sordo, R., Pailler, F., et al. 2023, A&A, 674, A26
  * Dantas et al. (2023) Dantas, M. L. L., Smiljanic, R., Boesso, R., et al. 2023, A&A, 669, A96
  * De Angeli et al. (2023) De Angeli, F., Weiler, M., Montegriffo, P., et al. 2023, A&A, 674, A2
  * de Sá-Freitas et al. (2023a) de Sá-Freitas, C., Fragkoudi, F., Gadotti, D. A., et al. 2023a, A&A, 671, A8
  * de Sá-Freitas et al. (2023b) de Sá-Freitas, C., Gadotti, D. A., Fragkoudi, F., et al. 2023b, arXiv e-prints, arXiv:2308.04482
  * Debattista et al. (2006) Debattista, V. P., Mayer, L., Carollo, C. M., et al. 2006, ApJ, 645, 209
  * Ellison et al. (2011) Ellison, S. L., Nair, P., Patton, D. R., et al. 2011, MNRAS, 416, 2182
  * Feuillet et al. (2018) Feuillet, D. K., Bovy, J., Holtzman, J., et al. 2018, MNRAS, 477, 2326
  * Fragkoudi et al. (2020) Fragkoudi, F., Grand, R. J. J., Pakmor, R., et al. 2020, MNRAS, 494, 5936
  * Fragkoudi et al. (2019) Fragkoudi, F., Katz, D., Trick, W., et al. 2019, MNRAS, 488, 3324
  * Freeman & Bland-Hawthorn (2002) Freeman, K. & Bland-Hawthorn, J. 2002, ARA&A, 40, 487
  * Fuhrmann (1998) Fuhrmann, K. 1998, A&A, 338, 161
  * Gaia Collaboration et al. (2023) Gaia Collaboration, Vallenari, A., Brown, A. G. A., et al. 2023, A&A, 674, A1
  * Grenon (1972) Grenon, M. 1972, in IAU Colloq. 17: Age des Etoiles, ed. G. Cayrel de Strobel & A. M. Delplace, 55
  * Guiglion et al. (2023) Guiglion, G., Nepal, S., Chiappini, C., et al. 2023, arXiv e-prints, arXiv:2306.05086
  * Halle et al. (2015) Halle, A., Di Matteo, P., Haywood, M., & Combes, F. 2015, A&A, 578, A58
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357
  * Helmi (2020) Helmi, A. 2020, ARA&A, 58, 205
  * Helmi et al. (2018) Helmi, A., Babusiaux, C., Koppelman, H. H., et al. 2018, Nature, 563, 85
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science and Engineering, 9, 90
  * Iles et al. (2023) Iles, E. J., Pettitt, A. R., Okamoto, T., & Kawata, D. 2023, arXiv e-prints, arXiv:2310.07011
  * Imig et al. (2023) Imig, J., Price, C., Holtzman, J. A., et al. 2023, ApJ, 954, 124
  * Isern (2019) Isern, J. 2019, ApJ, 878, L11
  * Kawata et al. (2018) Kawata, D., Baba, J., Ciucǎ, I., et al. 2018, MNRAS, 479, L108
  * Kawata et al. (2023) Kawata, D., Kawahara, H., Gouda, N., et al. 2023, arXiv e-prints, arXiv:2307.05666
  * Khoperskov et al. (2020a) Khoperskov, S., Di Matteo, P., Haywood, M., Gómez, A., & Snaith, O. N. 2020a, A&A, 638, A144
  * Khoperskov & Gerhard (2022) Khoperskov, S. & Gerhard, O. 2022, A&A, 663, A38
  * Khoperskov et al. (2020b) Khoperskov, S., Gerhard, O., Di Matteo, P., et al. 2020b, A&A, 634, L8
  * Kordopatis et al. (2015) Kordopatis, G., Binney, J., Gilmore, G., et al. 2015, MNRAS, 447, 3526
  * Leung & Bovy (2019) Leung, H. W. & Bovy, J. 2019, MNRAS, 489, 2079
  * Liang et al. (2023) Liang, X., Yoon, S.-J., Zhao, J., et al. 2023, The Astrophysical Journal, 956, 146
  * Lin et al. (2020) Lin, L., Li, C., Du, C., et al. 2020, MNRAS, 499, 1406
  * Lindegren et al. (2021) Lindegren, L., Bastian, U., Biermann, M., et al. 2021, A&A, 649, A4
  * Łokas (2020) Łokas, E. L. 2020, A&A, 634, A122
  * Matteucci (2001) Matteucci, F. 2001, The chemical evolution of the Galaxy, Vol. 253
  * Matteucci (2021) Matteucci, F. 2021, A&A Rev., 29, 5
  * McKinney (2010) McKinney, W. 2010, in Proceedings of the 9th Python in Science Conference, ed. Stéfan van der Walt & Jarrod Millman, 56 – 61
  * McMillan (2017) McMillan, P. J. 2017, MNRAS, 465, 76
  * Miglio et al. (2021) Miglio, A., Chiappini, C., Mackereth, J. T., et al. 2021, A&A, 645, A85
  * Minchev et al. (2013) Minchev, I., Chiappini, C., & Martig, M. 2013, A&A, 558, A9
  * Minchev et al. (2014) Minchev, I., Chiappini, C., Martig, M., et al. 2014, ApJ, 781, L20
  * Minchev & Famaey (2010) Minchev, I. & Famaey, B. 2010, ApJ, 722, 112
  * Minchev et al. (2011) Minchev, I., Famaey, B., Combes, F., et al. 2011, A&A, 527, A147
  * Monari et al. (2019) Monari, G., Famaey, B., Siebert, A., Wegg, C., & Gerhard, O. 2019, A&A, 626, A41
  * Montalbán et al. (2021) Montalbán, J., Mackereth, J. T., Miglio, A., et al. 2021, Nature Astronomy, 5, 640
  * Mor et al. (2019) Mor, R., Robin, A. C., Figueras, F., Roca-Fàbrega, S., & Luri, X. 2019, A&A, 624, L1
  * Pagel (1997) Pagel, B. E. J. 1997, Nucleosynthesis and Chemical Evolution of Galaxies
  * Palicio et al. (2023) Palicio, P. A., Recio-Blanco, A., Poggio, E., et al. 2023, A&A, 670, L7
  * Pérez-Villegas et al. (2017) Pérez-Villegas, A., Portail, M., Wegg, C., & Gerhard, O. 2017, ApJ, 840, L2
  * Pfenniger & Norman (1990) Pfenniger, D. & Norman, C. 1990, ApJ, 363, 391
  * Portail et al. (2017) Portail, M., Gerhard, O., Wegg, C., & Ness, M. 2017, MNRAS, 465, 1621
  * Queiroz et al. (2023) Queiroz, A. B. A., Anders, F., Chiappini, C., et al. 2023, A&A, 673, A155
  * Queiroz et al. (2018) Queiroz, A. B. A., Anders, F., Santiago, B. X., et al. 2018, MNRAS, 476, 2556
  * Queiroz et al. (2021) Queiroz, A. B. A., Chiappini, C., Perez-Villegas, A., et al. 2021, A&A, 656, A156
  * Recio-Blanco et al. (2023) Recio-Blanco, A., de Laverny, P., Palicio, P. A., et al. 2023, A&A, 674, A29
  * Riello et al. (2021) Riello, M., De Angeli, F., Evans, D. W., et al. 2021, A&A, 649, A3
  * Rocha-Pinto et al. (2000) Rocha-Pinto, H. J., Scalo, J., Maciel, W. J., & Flynn, C. 2000, ApJ, 531, L115
  * Sahlholdt et al. (2022) Sahlholdt, C. L., Feltzing, S., & Feuillet, D. K. 2022, MNRAS, 510, 4669
  * Schönrich et al. (2010) Schönrich, R., Binney, J., & Dehnen, W. 2010, MNRAS, 403, 1829
  * Sellwood (2014) Sellwood, J. A. 2014, Reviews of Modern Physics, 86, 1
  * Skrutskie et al. (2006) Skrutskie, M. F., Cutri, R. M., Stiening, R., et al. 2006, AJ, 131, 1163
  * Soderblom (2010) Soderblom, D. R. 2010, ARA&A, 48, 581
  * Steinmetz & Mueller (1994) Steinmetz, M. & Mueller, E. 1994, A&A, 281, L97
  * Sysoliatina & Just (2021) Sysoliatina, K. & Just, A. 2021, A&A, 647, A39
  * Taylor (2005) Taylor, M. B. 2005, in Astronomical Society of the Pacific Conference Series, Vol. 347, Astronomical Data Analysis Software and Systems XIV, ed. P. Shopbell, M. Britton, & R. Ebert, 29
  * Tepper-Garcia et al. (2021) Tepper-Garcia, T., Bland-Hawthorn, J., Vasiliev, E., et al. 2021, arXiv e-prints, arXiv:2111.05466
  * Trevisan et al. (2011) Trevisan, M., Barbuy, B., Eriksson, K., et al. 2011, A&A, 535, A42
  * Waskom (2021) Waskom, M. L. 2021, Journal of Open Source Software, 6, 3021
  * Wylie et al. (2022) Wylie, S. M., Clarke, J. P., & Gerhard, O. E. 2022, A&A, 659, A80
  * Xiang & Rix (2022) Xiang, M. & Rix, H.-W. 2022, Nature, 603, 599

##  Appendix A Details on input parameters for StarHorse and kinematic
calculations

The distances and stellar ages are computed with the StarHorse Bayesian
isochrone-fitting method (Queiroz et al. 2018, 2023) as noted in Sect. 2. As
inputs to StarHorse we use the spectroscopic parameters from G23, galactic
longitude (l) and latitude (b), photometric magnitudes G,
$\mathrm{\text{B}_{\text{p}}}$ $\&$ $\mathrm{\text{R}_{\text{p}}}$ and
parallaxes from _Gaia_ DR3 along with parallax corrections by Lindegren et al.
(2021). We also use the infra-red photometry (JHKs) from Two Micron All Sky
Survey (2MASS; Skrutskie et al. 2006).

To calculate positions and velocities in the galactocentric rest-frame and to
integrate the orbits of the stars, we use the 6D phase-space coordinates (sky
positions, parallaxes, proper motions and radial velocities) from _Gaia_ DR3
(Gaia Collaboration et al. 2023) along with the StarHorse distances. The
integration of orbits was done with Galpy (Bovy 2015), a python package for
Galactic dynamics calculations. We use Astropy (Astropy Collaboration et al.
2022) for coordinate and velocity transformations, assuming the Sun is located
at radius $\mathrm{R_{0}}$ = 8.2 kpc and the circular velocity of Local
Standard of Rest (LSR) as $\mathrm{V_{0}}$ = 233.1 km s-1 (Bovy 2015; McMillan
2017). The peculiar velocity of the Sun with respect to the LSR is
$\mathrm{(U,V,W)_{\odot}}=(11.1,\,12.24,\,7.25)$ km s-1(Schönrich et al.
2010). To run Galpy we adopt the MW potential of McMillan (2017) and perform
orbit integrations for a 3 Gyr period and save each orbit’s trajectory every 2
Myr.

The guiding radius we compute, given by $R_{g}=L_{Z}/V_{0}$, is independent of
the axisymmetric potential. Here, $L_{Z}$ denotes the star’s instantaneous
angular momentum, defined as $L_{Z}=R\cdot V_{\phi}$, where R is its
galactocentric distance, and $V_{\phi}$ denotes its azimuthal velocity in the
Galactic plane.

##  Appendix B Kinematics of the metal-rich stars

Figure 5: Velocity dispersions, radial ($\sigma_{V_{R}}$, dashed), azimuthal
($\sigma_{V_{\phi}}$, solid) and vertical ($\sigma_{V_{Z}}$, dotted), as a
function of stellar ages. Red circles represent the full sample while black
stars represent the metal-rich ([Fe/H]¿0.1) stars.

Here we compare the velocity dispersion of our entire sample (169 727 stars)
with that of stars with [Fe/H]¿0.1. We highlight two main points in Fig. 5,
namely: a) while $\sigma_{V_{\phi}}$ of the two samples are similar to each
other in the 4-9 Gyr range, metal-rich stars have systematically lower
$\sigma_{V_{R}}$ and $\sigma_{V_{Z}}$ (therefore systematically on cooler
orbits than the bulk of the disc stars); b) at around 3 Gyr metal rich stars
show a sudden increase in velocity dispersion.

We also notice that for ages above circa 9-10 Gyr the sample traced by the red
curves in Fig.5 shows a strong increase in velocity dispersion mainly due to
the dominant contribution of thick disc stars.

##  Appendix C Checking selection effects.

Figure 6: Kiel diagram and the CMD for the 256 metal-rich stars with 0.15 ¡
[Fe/H]¡ 0.3 and 2 ¡ Age (Gyr) ¡ 4. For illustration purpose, PARSEC stellar
isochrones, for ages 2, 3, 3.5, 4 & 5 Gyrs and [Fe/H]= 0.1, 0.2 and 0.3 dex,
are included. The CMD uses $JHKs$ absolute magnitudes from 2MASS which were
used as inputs for StarHorse to estimate ages. These stars are most probably
formed during the epoch of bar formation.

In order to clarify for the effects of any color/magnitude selection on the
stellar ages we show, in Fig. 6, Kiel diagram
($\mathrm{\textit{T}_{\text{eff}}}$ vs log(g)) and Color Magnitude Diagram
(CMD) for a sample of 256 metal-rich stars with 0.15 ¡ [Fe/H]¡ 0.3. PARSEC
stellar isochrones Bressan et al. (2012), which are also used by StarHorse,
have been overplotted. The plots show a very well-behaved distribution of the
metal-rich stars with no hard cuts on the color that hints at absence of young
metal-rich stars.

Given that our sample, obtained from Guiglion et al. (2023), is limited to
$\mathrm{\textit{T}_{\text{eff}}}$¡ 6300 K, and consists mostly of F, G, K
stellar types, one could assume that we are biased against younger stars (more
importantly younger SMR stars). Recently, Sahlholdt et al. (2022), using a
comparable data set of $\sim$ 180 000 MSTO and SGB stars from the GALAH DR3
(Buder et al. 2021), study the age-metallicity distribution. Their sample
covers a wider 8000 ¡ $\mathrm{\textit{T}_{\text{eff}}}$ ¡ 4000 range
including hotter stars not present in our RVS-CNN sample (see their Fig. 1).
Their AMR (see their Fig. 2) also shows that an absence of SMR stars younger
than 2-3 Gyrs is real and not an effect of selection bias. This further
supports our result on a dearth of SMR stars younger than 3 Gyrs ago as
discussed in Sect. 3.2. Also see Appendix D.

##  Appendix D Validation with external catalog.

Here, we perform the validation of our finding, of the steep decline of upper
[Fe/H] envelope starting $\sim$3 Gyrs ago as discussed in Sect. 3.2, using an
external catalog. We use the Red-giant sample in the Kepler field from the
Miglio et al. (2021) paper, reanalyzed with APOGEE DR17 parameters. After
recommended quality cuts and selecting stars with age ¡ 12 Gyr and
[Fe/H]¿-0.7, we obtained a sample of 2 614 stars with mean age uncertainty of
25%. Figure 7 shows the time evolution of the upper envelope of the [Fe/H] for
the asteroseismic sample. Similar to Fig. 3, the plot shows that our Galaxy
attained high levels of [Fe/H] enrichment quite early with a slow decline
until $\sim$ 3 Gyrs ago and then a steep decline after that. This confirms our
finding of the dearth of SMR stars younger than 3 Gyrs old in the solar
vicinity. In Fig 8 we show that the metal rich stars ([Fe/H]¿0.1) also show a
bimodality in the guiding radius, similar to the finding with with the RVS-CNN
sample as discussed in Sect. 3.1.

Figure 7: Upper envelope of the age-metallicity relation traced for the stars
with asteroseismic ages from Miglio et al. (2021) and [Fe/H] from the APOGEE
DR17 catalog (Abdurro’uf et al. 2022). The triangles represent the median
[Fe/H] of 5 (red) / 10 (black) most-metal-rich stars in each age bin and the
error-bars show the scatter in [Fe/H]. Smaller number of stars per bin, for
median and scatter, are used due to smaller number of stars in the Miglio et
al. (2021) catalog. Figure 8: Distribution of the guiding radius
($\mathrm{R_{g}}$) for the sample of metal-rich ([Fe/H]¿0.1) red-giant stars
with asteroseismic ages from the Miglio et al. (2021) catalog.

##  Appendix E Guiding Radii vs Ages

In Fig. 9, we present the stellar ages as a function of guiding radii, similar
to left plot of Fig. 3, in [Fe/H] bins similar to Fig. 2. The [Fe/H] values
for each bins increases as we move from panel a) to panel i). The plot shows
large number of stars, with a significant presence of stars younger than
$\sim$3 Gyrs, for the sub-solar to solar metallicity bins while the
$\mathrm{R_{g}}$ bimodality is not clear. As we move to [Fe/H]¿0.1, the
$\mathrm{R_{g}}$ bimodality becomes stronger with increasing [Fe/H]. The
figure shows an absence metal rich ([Fe/H]¿0.1) stars younger than 2.7 Gyrs,
marked by the dashed line. We also see that the most metal-rich stars
([Fe/H]¿0.3) are mostly old (11-6 Gyrs).

Figure 9: Stellar ages as a function of guiding radii, similar to left plot of
Fig. 3, in [Fe/H] bins similar to Fig. 2. A dashed black line at 2.7 Gyr is
also shown. The colors represent number of stars per bin. A Kernel Density
Estimate (KDE) has also been overplotted. As the [Fe/H] increases we find the
bimodality in guiding radius appear more clearly. At solar and sub-solar
metallicities we find a lot of young stars while for [Fe/H]¿0.1 we see a near
absence of stars younger than 2.7 Gyr old.

##  Appendix F Metal-rich stars as traces of bar activity

In Sect. 3.1 we discussed the [Fe/H] dependence of $\mathrm{R_{g}}$ bimodality
revealed by the metal-rich stars. Here we test further to show how the metal-
rich stars can trace bar activity.

Substructures and ridges found in the SNd have been well studied and
attributed to the bar (for detailed discussions see Kawata et al. (2018);
Fragkoudi et al. (2019); Chen et al. (2022)). Similarly, the Hercules Moving
Group, which shows higher metallicity compared to the disc, has been proposed
to consist of stars undergoing radial migration (Liang et al. 2023) and
corotation of the bar has been suggested as the physical mechanism for its
formation (e.g. Pérez-Villegas et al. 2017; Monari et al. 2019; Khoperskov &
Gerhard 2022; Liang et al. 2023).

In Fig. 10, we show our stars in the R vs $V_{\phi}$ plane to reveal the
substructures in the phase-space. The density distribution, in top panels,
clearly shows ridges in the phase space for the whole sample (thin disc stars
with $Z_{max}$¡1 kpc). These ridges are not visible for the metal-poor stars
(mid panel) while they are clearly seen for the metal-rich stars. Furthermore,
in the bottom panels we show how the undulations in $V_{R}$ are traced by our
sample, primarily due to precise distances thanks to Gaia. These undulations
are very clear for the full sample and for the metal-rich set, however the
signature is noisy for the metal-poor sample.

In Fig. 10, in middle and right panels, a black dashed line is drawn to show
the location of the top sharp ridge seen in the full sample (top-left). We
observe an interesting feature showing a near absence of metal-rich stars
beyond this line while metal-poor stars are present. Khoperskov & Gerhard
(2022) using high-resolution N-body simulation, attribute this sharp decrease
in mean [Fe/H] to the effect of OLR of the bar which limits the migration of
the metal-rich stars, coming from inner Galaxy, beyond OLR.

These tests show that the metal-rich stars are indeed excellent tracers for
bar activity.

Figure 10: Top row: 2D density distribution showing the azimuthal velocity
($V_{\phi}$) as a function of galactocentric radius (R) for the disc stars
(i.e. $Z_{max}<1$ kpc) to reveal the ridge-like features. Left panel shows all
stars, mid panel shows metal-poor stars ([Fe/H]¡$-$0.2) and right panel the
metal-rich stars ([Fe/H]¿0.1). Bottom row: Same plots as above shown to reveal
the undulations when colored by the $V_{R}$. Both ridges and undulations are
clearly seen for the full sample and the metal-rich stars while these features
are not so clear for the metal-poor stars. In middle and right panels a dotted
black line is drawn to represent the top sharp ridge seen in the top-left
density plot.