The Unmixed Debris of Gaia-Sausage/Enceladus in the Form of a Pair of Halo Stellar Overdensities

In the first billion years after its formation, the Galaxy underwent several mergers with dwarf satellites of various masses. The debris of Gaia-Sausage/Enceladus (GSE), the galaxy responsible for the last significant merger of the Milky Way, dominates the inner halo and has been suggested to be the progenitor of both the Hercules-Aquila Cloud (HAC) and Virgo Overdensity (VOD). We combine SEGUE, APOGEE, Gaia, and StarHorse distances to characterize the chemodynamical properties and verify the link between HAC, VOD, and GSE. We find that the orbital eccentricity distributions of the stellar overdensities and GSE are comparable. We also find that they have similar, strongly peaked, metallicity distribution functions, reinforcing the hypothesis of common origin. Furthermore, we show that HAC and VOD are indistinguishable from the prototypical GSE population within all chemical-abundance spaces analyzed. All these evidences combined provide a clear demonstration that the GSE merger is the main progenitor of the stellar populations found within these halo overdensities.


INTRODUCTION
The concordance cosmological model of galaxy formation predicts that the halos of massive galaxies, such as the Milky Way, are formed via the accretion of many low-mass satellites (Springel et al. 2006). While the star-by-star characterization of ancient accretion events in other galaxies is still unfeasible, wide-area photometric and spectroscopy surveys have shown the complexity of the Milky Way's stellar halo (or simply "halo"; Ivezić et al. 2012;Belokurov 2013) which includes, e.g., several stellar overdensities (Newberg & Carlin 2016). These overdensities are identified as stellar count excesses in given Galactic regions when compared to other homogeneous/smooth halo fields.
The stellar halo overdensities can form via the buildup of tidal debris, from one or more satellite galaxies, whose stars have highly eccentric orbits that accumulate at the apocenter (Newberg & Carlin 2016;Li et al. 2016;Deason et al. 2018), a scenario that has been explored with pure N -body models (Johnston et al. 2008;Helmi et al. 2011;Naidu et al. 2021). Moreover, the accretion of a single dwarf galaxy can create more than a single stellar overdensity, challenging the association between distinct substructures and their progenitors (e.g., Johnston et al. 2012). This hypothesis has gained recent attention with the confirmation of the last significant merger (mass ratio 1:4) experienced by the Milky Way with a dwarf galaxy named Gaia-Sausage/Enceladus (GSE; Belokurov et al. 2018;Helmi et al. 2018). The GSE merger happened ∼10 Gyr ago (Gallart et al. 2019;Bonaca et al. 2020; Montalbán et al. 2021) and is likely responsible for many features observed in our Galaxy Di Matteo et al. 2019;Iorio & Belokurov 2019;Belokurov et al. 2020 Currently, there are two overdensities that have been tentatively linked to GSE, the Hercules-Aquila Cloud (HAC) and the Virgo Overdensity (VOD). HAC was discovered as a diffuse substructure of main-sequence turnoff stars located toward the Galactic center (Belokurov et al. 2007). It extends in heliocentric distances (d ) from 10 to 20 kpc, from 25 • < l < 60 • and −40 • < b < 40 • (Belokurov et al. 2007;Watkins et al. 2009;Sesar et al. 2010;Simion et al. 2014). Watkins et al. (2009) and Sesar et al. (2010) VOD was first identified with RR Lyrae and mainsequence turnoff stars as an excess in the stellar halo (Vivas et al. 2001;Newberg et al. 2002). Later, these were associated as part of a larger overdensity located in the direction of the Virgo constellation (Jurić et al. 2008). Its main stellar excess is located between 270 • < l < 330 • and 50 • < b < 75 • , and possibly extends over a larger area on the sky (∼3000 deg 2 ; Jurić et al. 2008;Bonaca et al. 2012). Its [Fe/H] varies from −1.1 to −1.95 depending on the stellar sample (Vivas et al. 2008;Carlin et al. 2012;Naidu et al. 2021), and the stars associated with VOD are typically on high-eccentricity orbits ( e ∼ 0.8; Carlin et al. 2012;Simion et al. 2019). Simion et al. (2019) explored the association between the GSE merger and both VOD and HAC using astrometric data from the Gaia space mission (Gaia Collaboration et al. 2016). They showed that these overdensities share similar dynamical properties, such as the orbital energy, apocentric distance, and eccentricity, with GSE stars and interpreted it as an evidence that HAC and VOD are unmixed debris of GSE. This is also suggested by the apparent orbital pileup of local halo stars in the regions associated with these overdensities (Balbinot & Helmi 2021) and with pure N -body models mimicking the GSE merger (Naidu et al. 2021).
Due to their relative large distances from the Sun, it is difficult to obtain, for a large number of stars, accurate abundances and reliable distance measurements for members of these substructures. We overcome this challenge by taking advantage of the StarHorse code (Santiago et al. 2016;Queiroz et al. 2018), which is as a Bayesian isochrone-fitting method able to derive distances, extinctions, and other stellar parameters based on observed spectroscopic, photometric, and astrometric data. Crucially, these precise StarHorse distance estimates allows to confidently select members of these overdensities in order to conduct a chemical comparison between HAC/VOD and nearby GSE stars. This exercise is necessary to confirm (or reject) the link between the overdensities and the main accretion event of our Galaxy. Our goal is to verify whether i ) they have solely an accreted origin, ii ) they have similar chemo-dynamical properties among them, and iii ) they are (in)distinguishable from GSE.
In this spirit, we carry out an in-depth chemodynamical characterization of HAC and VOD with data from Apache Point Observatory Galactic Evolution Experiment (APOGEE; Majewski et al. 2017) and Sloan Extension for Galactic Understanding and Exploration (SEGUE; Yanny et al. 2009). With the new data, we are now able to demonstrate that the majority of the stars from the stellar overdensities are consistent with the prototypical population of GSE defined within 5 kpc around the Sun. We keep only those stars from SEGUE (Rockosi et al. 2022) with S/N > 20 pixel −1 . We limit our SEGUE sample to [Fe/H] < −0.5 to decrease contamination from disk stars. We also select stars within 4500 < T eff /K < 6500, which ensures the stellar parameters are accurately determined and avoids large uncertainties. The abundance data from SEGUE are less precise (uncertainties in [Fe/H] ≥ 0.1), but with the advantage of observing toward fainter-magnitude stars, providing a larger sample for each overdensity, which allows a robust determination of the orbital eccentricity and metallicity distributions. To obtain reliable abundances for a larger set of elements, we expand our analysis using APOGEE data.
We combine (1.5" search radius) the aforementioned samples with Gaia Early Data Release 3 (Gaia Collaboration et al. 2020) to obtain the absolute proper motions and their uncertainties. To ensure good astrometric solutions, we impose renormalized unit weight errors within the recommended range (RUWE ≤ 1.4; Lindegren et al. 2020). Besides that, we remove stars with parallax < −5 and with low fidelity (fidelity v2 < 0.5; see Rybizki et al. 2022). The radial velocities are from the spectroscopic surveys, and the distances were obtained with the StarHorse code, specifically from the new releases where several large-scale spectroscopic surveys are matched to Gaia EDR3 data. For APOGEE DR17, distance estimates are publicly available 5 and, for SEGUE, these will be made available alongside a forthcoming paper (Queiroz et al., in preparation). We selected only stars with fractional (Gaussian) uncertainties of their nominal distance values <20%.
We distinguish HAC into two substructures HAC-South (HAC-S) and HAC-North (HAC-N). Under the assumption that both regions have the same origin (Iorio & Belokurov 2019;Simion et al. 2019;Naidu et al. 2021), they should have similar chemodynamical properties. We select stars from HAC-S and HAC-N with the following criteria: 30 • < l < 60 • , −45 • < b < −20 • and 20 • < b < 45 • , respectively. For VOD, our selection is based on the following Galactic coordinate cuts: 270 • < l < 330 • and 50 • < b < 75 • . For both overdensities, we consider d between 10 and 20 kpc, with a mean fractional uncertainty of ∼13% in this range. These delineated criteria were applied to both SEGUE and APOGEE samples. Figure 1 shows the resulting spatial distribution of stars selected from HAC and VOD in SEGUE.

Orbital Parameters
To obtain the orbital parameters, we employ the axisymmetric Galactic potential of McMillan (2017). The adopted distance from the Sun to the Galactic center is 8.2 kpc (Bland-Hawthorn & Gerhard 2016), the circular velocity at this position is v circ = 232.8 km s −1 , and the peculiar motion of the Sun with respect to the local circular orbit is (U, V, W ) = (11.10, 12.24, 7.25) km s −1 (Schönrich et al. 2010). We integrated the orbits forward in time for 10 Gyr. For each star, we construct a set of 100 initial conditions using a Monte Carlo technique taking into account the observational uncertainties in d , proper motions, and the line-of-sight velocities. The final dynamical parameters are taken as the medians of the derived distributions and associated uncertainties are the 16 th and 84 th percentiles.
In order to compare whether the stellar overdensities share similar dynamical properties with GSE, we use a chemodynamical criteria (Limberg et al. 2022; see their equation 1) 6 to identify nearby-within 5 kpc from the Sun-GSE stars, which is designed to yield maximum purity. The final SEGUE/APOGEE samples of HAC-S, HAC-N, and VOD contain 231/11, 516/21, and 556/22 stars, respectively. 6 We note that, since SEGUE does not provide abundances other than [Fe/H] and [α/Fe], we do not apply the chemical portion of these author's criteria (which demand Mg, Mn, and Al).

. Eccentricity and Metallicity Distributions
We take advantage of the large number of the HAC-S, HAC-N, and VOD stars found in SEGUE to characterize the orbital eccentricity and metallicity distributions of these stellar overdensities with the purpose of testing their connection with GSE. We explore their distributions of orbital eccentricity in Figure 2. The first panel shows the full sample, and overall, the overdensities have a similar eccentricity distributions; they are dominated (∼65%) by stars with e > 0.7, similar to GSE Myeong et al. 2018;Limberg et al. 2021), but have an extended tail toward lower eccentricity values.
We also show the eccentricity distribution at different metallicity intervals (second to fifth panels of Figure 2). The overdensities have similar distributions across all the metallicity intervals but differ from the GSE. For instance, in the very metal-poor regime (VMP; [Fe/H] < −2.0), the overdensities have roughly equal amounts of stars in each eccentricity bin, whereas the GSE continues to be dominated by stars with e > 0.7. Toward higher metallicities (−1.5 < [Fe/H] < −1.0), where the peak of the GSE metallicity distribution function (MDF) is located (see Limberg et al. 2022 and references therein), the overdensities are dominated by stars with e > 0.7 and, differently from GSE, show an extended tail toward lower eccentricities.
In addition, Figure 2 shows that the contribution of low-eccentricity stars is higher toward the lowmetallicity regime ([Fe/H] < −1.5). The absence of the low-eccentricity tail in the GSE sample could be a bias due to our selection criterion, which favors stars with vertical component of the angular momentum (L z ) ∼ 0. On the other hand, simulations mimicking GSElike mergers (Naidu et al. 2021;Amarante et al. 2022) show that the vast majority of the accreted stars end up with e > 0.7. Indeed, models that take into account star formation and a self-consistent chemical en-richment show that the accreted stars on relatively lower eccentricities are more metal-poor compared to those on higher eccentricities (see figure 6 in Amarante et al. 2022). Only a small fraction (∼1−8%) of the accreted stars end up on less eccentric orbits (e < 0.7) and these are typically more metal-poor (by ∼0.2−0.7 dex) than those with e > 0.7. Therefore, part of the metal-poor stars on low-eccentricity orbits found in SEGUE could be the metal-poor tail of GSE, as expected by the aforementioned idealized GSE-like merger models. We note that the metal-poor stars in our HAC/VOD sample have −1500 L z 1500 kpc km/s and have roughly the same fraction of prograde and retrograde orbits. Thus, other mildly eccentric (e ∼ 0.5) exclusively retrograde (L z −1500 kpc km/s) halo structures can be discarded as responsible for the low-metallicity and loweccentricity component identified here.
The presence of the low-eccentricity stars in both HAC-S, HAC-N, and VOD may also be interpreted as the existence of other merger(s) event(s) that contributed less significantly to the formation of these overdensities. We verified that the upper limit of contamination by stars from the Sagittarius stream in both overdensities is <5% following the criteria of Naidu et al. (2020) (see also Peñarrubia & Petersen 2021). Another potential contributor to the low-eccentricity regime might be LMS-1/Wukong (Yuan et al. 2020;Naidu et al. 2020;Horta et al. 2022), which contains mostly VMP stars with e < 0.7. This structure is also composed of stars on polar orbits (J z > J r ) that spatially overlap with VOD and HAC (Yuan et al. 2020). We note, however, that LMS-1 members are exclusively on prograde orbits, which is not the case for the VMP low-eccentricity stars in our HAC/VOD sample.
Finally, in situ stellar populations could contribute to the stellar overdensities, such as the Splash (Belokurov et al. 2020) and Aurora (Belokurov & Kravtsov 2022). The Splash consists of metal-rich ([Fe/H] > −0.7) stars on mildly eccentric orbits (e > 0.5), while the ma- jority of HAC and VOD low-eccentricity stars have [Fe/H] < −1.5. Aurora is mainly characterized by its thick-disk-like chemistry, i.e., high [α/Fe] at the same [Fe/H] in comparison to the Galactic halo/GSE. Moreover, Aurora should be identifiable at metallicities of −1.5 < [Fe/H] < −1.0, which is not the case for the loweccentricity populations identified within HAC/VOD. Therefore, given the present data, we conclude that is is unlikely that an in situ contribution could be the major non-GSE contaminant to the stellar overdensities. Figure 3 shows the MDF for the full sample (left panel) and different eccentricity intervals (middle and right panels) for the SEGUE stars from the overdensities and GSE. In the left panel, the MDF of GSE is similar to those from the different overdensities, with a slight difference toward the VMP regime. An equivalent result has been obtained by Naidu et al. (2021), with metallicities from the H3 survey.
In the middle panel, we show the MDFs for stars with e > 0.7, which is the range of the eccentricity characteristic of GSE. We estimate the medians and median absolute deviations of the distributions by bootstrap sampling with 10 4 iterations, accounting for the uncertainties in metallicity. The GSE median [Fe/H] value in the SEGUE sample is −1.39 −0.01 +0.01 dex with median absolute deviation of 0.20 dex, which is in agreement with the median [Fe/H] of HAC-S, HAC-N, and VOD; −1.51 −0.06 +0.06 , −1.47 −0.05 +0.05 , and −1.35 −0.03 +0.03 dex, respectively. Differences in the median [Fe/H] between all overdensities and GSE are smaller than the typical uncertainties in SEGUE (≥0.1 dex). We also note that these listed median [Fe/H] values are well within literature determinations for GSE (see discussion in Limberg et al. 2022).
In the right panel of Figure 3, we show the MDFs for e < 0.7. The most noticeable feature in this loweccentricity regime is the excess of VMP stars in both HAC and VOD in comparison to their high-eccentricity counterparts (middle panel). This behavior is the same as previously identified in Figure 2. The remaining stars from our local GSE selection (<1%) are also mostly VMP, which might indicate either contamination from other accreted substructures (Helmi 2020) or a lowmetallicity component of GSE itself, although we reinforce that the lack of low-eccentric GSE stars might be a bias due to our selection criteria.

Elemental Abundances and Orbital Parameters
We took advantage of the SEGUE sample to show that the bulk of the stars from the stellar overdensities have similar eccentricity and metallicity distributions, which are also compatible with GSE, suggesting that they share a common origin. We now explore in detail the chemical-abundance patterns of the overdensities with the available data from APOGEE. Figure 4 shows the abundance of stars from HAC-S, HAC-N, VOD, and GSE  In all panels, the yellow lines are isodensity contours associated with our GSE sample and the background 2D histograms show stars from APOGEE. All the overdensities clearly overlap with the contour of GSE and occupy the region of the chemical abundance planes dominated by accreted populations. This indicates that HAC-S/N and VOD share chemical properties compatible with GSE, favoring a scenario of a common origin for these substructures. However, we note that the origin of the VMP stars on low-eccentricity orbits identified in Section 3.1 could not be explored with APOGEE data given that it does not reach such low-metallicity regime.
To further analyze this chemodynamical compatibility, we select stars from HAC and VOD that have the characteristics of ancient accreted populations in the  Hawkins et al. 2015;Das et al. 2020) or have low-eccentricity orbits (e < 0.7). The yellow contours indicate the chemical patterns of the GSE stars. The background 2D histogram shows the distribution of the APOGEE stars.
[Mg/Mn]-[Al/Fe] space (Hawkins et al. 2015;Das et al. 2020) and eccentricity consistent with GSE members (e > 0.7). The stars that do not satisfy these criteria (region of [Mg/Mn]-[Al/Fe] diagram delimited by red lines; Limberg et al. 2022) are represented by the empty symbols in Figure 4, and we note that only a small fraction (∼12%) occupy the in situ locus. This exercise makes it clear that the bulk of stars of the stellar overdensities studied in this work occupy the same region of GSE population in all chemical-abundance planes.
In Figure 5, we show orbital parameters of HAC-S, HAC-N, VOD, and GSE. In the left panel, we see that the stellar overdensities occupy a well-defined region, which is dominated by GSE stars, in the total orbital en-ergy (E) versus L z plane. The right panel (e-L z ) shows that the majority of HAC/VOD stars in this sample are on highly eccentric orbits and low average L z , which are dynamical properties of GSE. In general, all the aforementioned chemodynamical properties suggest that the halo overdensities studied in this work share a common origin with GSE.
The stellar overdensities of our study are dynamically compatible with the GSE, which dominates the stellar content of the inner halo Naidu et al. 2020). This combination of orbital parameters was also identified in other works (Simion et al. 2019;Balbinot & Helmi 2021) and is in contradiction with the hypothesis that these structures are dominated by stars on polar orbits or connected in a polar ring as speculated by Jurić et al. (2008).

CONCLUSIONS AND FINAL REMARKS
In this Letter, we combined SEGUE, APOGEE, and Gaia data, together with StarHorse distances, to investigate the origin of both HAC and VOD and their potential link with GSE. With this goal, we performed a detailed chemodynamical analysis of these overdensities and compared them with GSE. Our main results are summarized below.
• The stellar overdensities studied in this work show similar eccentricity distributions between them. VOD and HAC are composed mostly of stars with e > 0.7 (68% and 65%, respectively), which is characteristic of GSE.
• • We identified that HAC-S, HAC-N, and VOD exhibit contributions of metal-poor stars ([Fe/H] −1.5) on low-eccentricity orbits (e ∼ 0.5). We speculate some scenarios for the origin of these stars as either a less eccentric and metal-poor portion of GSE or the contribution of other accretion event(s). We do not rule out the contribution/contamination of in situ stars in the stellar overdensities, although this is difficult to reconcile with their low [Fe/H] values.
Our main result is that the majority of stars from HAC and VOD have chemodynamical properties characteristic of accreted populations, which are similar between them and indistinguishable from GSE. This serves as a constraint to models that simulate this accretion event and should reproduce stellar overdensities such as HAC and VOD in the inner Galactic halo. Furthermore, additional high-resolution spectroscopic observations of HAC and VOD members with different eccentricities will allow us to understand the origin of the group of VMP stars on low-eccentricity orbits.