key: cord-0547751-euro7pb0
authors: Horta, Danny; Schiavon, Ricardo P.; Mackereth, J. Ted; Pfeffer, Joel; Mason, Andrew C.; Kisku, Shobhit; Fragkoudi, Francesca; Prieto, Carlos Allende; Cunha, Katia; Hasselquist, Sten; Holtzman, Jon; Majewski, Steven R.; Nataf, David; O'Connell, Robert W.; Schultheis, Mathias; Smith, Verne V.
title: Evidence from APOGEE for the presence of a major building block of the halo buried in the inner Galaxy
date: 2020-07-20
journal: nan
DOI: nan
sha: 355d720537bd58e5d35e3b186e602f4288f33870
doc_id: 547751
cord_uid: euro7pb0

We report evidence from APOGEE for the presence of a new metal-poor stellar structure located within $sim$4~kpc of the Galactic centre. Characterised by a chemical composition resembling those of low mass satellites of the Milky Way, this new inner Galaxy structure (IGS) seems to be chemically and dynamically detached from more metal-rich populations in the inner Galaxy. We conjecture that this structure is associated with an accretion event that likely occurred in the early life of the Milky Way. Comparing the mean elemental abundances of this structure with predictions from cosmological numerical simulations, we estimate that the progenitor system had a stellar mass of $sim5times10^8M_odot$, or approximately twice the mass of the recently discovered Gaia-Enceladus/Sausage system. We find that the accreted:{it in situ} ratio within our metal-poor ([Fe/H]$<$--0.8) bulge sample is somewhere between 1:3 and 1:2, confirming predictions of cosmological numerical simulations by various groups.

sought to characterise halo stellar populations, associating them with in situ formation or an accretion origin. Detection of substructure in phase space has worked very well for the identification of ongoing and/or recent accretion events (e.g., Ibata et al. 1994; Helmi et al. 1999; Majewski et al. 2003; Belokurov et al. 2006) . However, phase mixing makes the identification of early accretion activity that occurred during the first few Gyr of the Milky Way's life more difficult, requiring additional information, usually in the form of detailed chemistry (e.g., Nissen & Schuster 2010; Hayes et al. 2018; Mackereth et al. 2019 ).

The combination of data from the Gaia satellite (Gaia Collaboration et al. 2018 ) and chemistry from massive high resolution spectroscopic surveys (e.g., Gilmore et al. 2012; Majewski et al. 2017; Martell et al. 2017) is transforming this field. In the past few years, several distinct satellite accretions have been suggested by various groups, including the Gaia-Enceladus/Sausage system Helmi et al. 2018; Belokurov et al. 2018; Mackereth et al. 2019) , Sequoia (Myeong et al. 2019) , the Kraken (Kruijssen et al. 2019b) , and Thamnos 1 and 2 . More recently, a number of additional structures were identified by Naidu et al. (2020) on the basis of data collected as part of the H3 survey (Conroy et al. 2019) . While the reality of some of those events still needs to be fully established (e.g., Jean-Baptiste et al. 2017; Koppelman et al. 2020) , it certainly is clear that we are just scratching the surface and the field is ripe for exciting new findings in the near future.

The central few kpc of the Galactic halo are obviously extremely important when it comes to retelling the early accretion history of the Milky Way and discerning the contribution of in situ formation to the stellar halo mass. Assuming a single power law density profile with exponent α = −2.96 (Iorio et al. 2018 , Horta et al. 2020 with spherical symmetry, roughly 50% of the halo stellar mass (out to R GC ∼30 kpc) is located within ∼3 kpc of the Galactic centre. Moreover the inner few kpc of the Galaxy is the region one would expect to host most of the early in situ halo star formation, including the oldest stars in the Galaxy (e.g., Tumlinson 2010; Savino et al. 2020) , as well as the remnants of early accretion events driven there by dynamical friction (e.g., Tremaine et al. 1975; Pfeffer et al. 2020) . Observational access to inner halo populations is however quite difficult. They inhabit a region of the Milky Way that by convention is referred to as the Galactic bulge, whose stellar population content consists of a superposition of structures that are far more densely populated than the halo itself, including the thick and thin disks, and the bar (e.g., Ness et al. 2013a,b; Rich 2013; Nataf 2017 ; Barbuy et al. 2018) . Moreover, observational access is made difficult by severe extinction towards the inner several degrees from the Galactic centre.

Since the groundbreaking work by Rich (1988) , a truly herculean effort has been invested by several groups towards mapping the stellar population content of the Galactic bulge, mostly based on optical spectroscopy (e.g., McWilliam & Rich 1994; Ibata et al. 1994; Zoccali et al. 2006; Fulbright et al. 2007; Howard et al. 2008; Gonzalez et al. 2011; Rich et al. 2012; García Pérez et al. 2013; Rojas-Arriagada et al. 2014; Ryde et al. 2016; Rojas-Arriagada et al. 2017) 1 . However, until very recently, the detailed chemistry of the metalpoor population inhabiting the inner few kpc of the Galactic centre has been poorly characterized, with most studies being based on small samples with metallicities below [Fe/H] ∼ −1. With the advent of the Apache Point Observatory Galactic Evolution Experiment (Majewski et al. 2017) , detailed chemistry (combined with precision radial velocities) has recently been obtained for a statistically significant sample of metal-poor stars within a few kpc of the Galactic centre. A number of studies resulted from that still growing database. Nidever et al. (2012) identified cold velocity peaks in the radial velocity distribution of stars within ∼ 20 • of the Galactic centre, which may be associated with resonant bar orbits (Molloy et al. 2015) . García Pérez et al. (2013) identified some of the most metal-poor stars known in the Galactic bulge; Zasowski et al. (2016) showed that the kinematics of metal-poor bulge stellar populations is dynamically hotter than that of their metal-rich counterparts; Schultheis et al. (2017) found evidence for the presence of a young component among the bulge stellar populations, in agreement with previous work; Schiavon et al. (2017) discovered a large number of field stars likely resulting from the destruction of an early population of globular clusters; García Pérez et al. (2018) studied the metallicity distribution function (MDF) of bulge populations; Zasowski et al. (2019) scrutinised the detailed chemistry of bulge populations, showing that they seem to be consistent with a single evolutionary track, across a wide metallicity range; and Queiroz et al. (2020, in prep.) made a detailed characterization of bulge stellar populations on the basis of their kinematics and chemistry.

In this paper we report evidence based on APOGEE spectroscopy, Gaia DR2 distances, and EAGLE numerical simulations, for the presence of a metal-poor stellar population within ∼4 kpc of the Galactic centre, that is both chemically and dynamically distinct from its co-spatial, more metal-rich counterparts. We argue that this population is likely the remnant of a massive accretion event that probably occurred in the early stages of the Milky Way formation. Due to its location at the heart of the Galaxy, we name this system the Inner Galaxy Structure (IGS). In Section 2 we briefly describe the data upon which this work is based, and our sample selection criteria. In Section 3 we discuss the chemical criteria to distinguish accreted from in situ populations. In Section 4 we discuss the identification of substructure in integrals of motion (IoM) space, arguing that the newly identified structure is dynamically detached from other components of the inner Galaxy. Section 5 presents a discussion of the chemical properties of the new structure, where we provide evidence that it is chemically detached from the other populations inhabiting the inner Galaxy. In Section 6 we discuss the origins of the metal-poor stellar populations providing an assessment of the contribution of accreted and in situ populations to the stellar mass budget of the inner halo. We also estimate physical properties of the putative satellite progenitor of the IGS, such as mass and star formation history. Our conclusions are summarised in Section 7.

This paper is based on a combination of data from the SDSS/APOGEE survey (Blanton et al. 2017; Majewski et al. 2017) , data release 16 (DR16, Ahumada et al. 2019) , with 6D phase information and orbital parameters inferred from Gaia DR2 distances and proper motions (Gaia Collaboration et al. 2018 ) and APOGEE radial velocities (Nidever et al. 2015) . We make use of the distances from Leung & Bovy (2019b) for the APOGEE DR16 data, generated using the astroNN python package (for a full description see Leung & Bovy 2019a ). These distances are determined using a previously trained astroNN neural-network, which predicts stellar luminosity from spectra using a training set comprising of stars with APOGEE spectra and Gaia DR2 parallax measurements (Gaia Collaboration et al. 2018) . The model is able to simultaneously predict distances and account for the parallax offset present in Gaia-DR2, producing high precision, accurate distance estimates for APOGEE stars, which match well with external catalogues and standard candles. Orbital parameters were calculated using the publicly available code galpy 2 (Bovy 2015; Mackereth & Bovy 2018) , and employing a McMillan (2017) potential. We checked the quality of the Gaia astrometry for the whole sample finding that, save for a very small number of outliers, the renormalized unit weight error (RUWE) for our bulge sample is within the safe recommended region (RUWE<1.44) 3 , which speaks for the reliability of our orbital parameters.

APOGEE data are based on observations collected with twin high resolution multi-fiber spectrographs (Wilson et al. 2019) attached to the 2.5 m Sloan telescope at Apache Point Observatory (Gunn et al. 2006 ) and the du Pont 2.5 m telescope at Las Campanas Observatory (Bowen & Vaughan 1973) . Elemental abundances are derived from automatic analysis of stellar spectra using the ASPCAP abundance pipeline (García . The spectra themselves were reduced by adoption of a customized pipeline (Nidever et al. 2015) . For a detailed description of the APOGEE data, see Jönsson et al. (2020) for DR16 and Holtzman et al. (2015, 2018) and Jönsson et al. (2018) for previous data releases. For details on targeting procedures, see Zasowski et al. (2017) .

The parent sample upon which our work is based is defined as follows. Stars were selected from the DR16 catalog that match the following criteria: 4000 < T eff < 6000 K, 1 < log g < 3, [Fe/H]> −1.7, S/N>70, d err /d < 0.2 (where d 2 https://github.com/jobovy/galpy 3 See discussion in https://www.cosmos.esa.int/web/gaia/dr2known-issues [Al/Fe] plane. The 2D histogram shows the entire parent sample, and the red dots highlight stars located at R GC < 4 kpc. The black solid line defines our criterion to distinguish in situ from accreted populations. The dotted line further splits the latter group between low-and high-α stars. See details in Section 2. and d err are distance and its uncertainty). The S/N and stellar parameter criteria are designed to maximise the quality of the elemental abundances, minimising systematic effects at low and high T eff , minimising/eliminating contamination by AGB and nearby dwarf stars at low and high log g, respectively, and stars with very low metallicity, for which abundances of elements such as Mn and N become uncertain. We also remove stars for which ASPCAP did not provide reliable abundances for Fe, Mg, Al, and Mn, and any stars with STAR_BAD flags set. Finally, 6,022 globular cluster stars, identified following the procedure described by Horta et al. (2020) , were also removed from consideration. The resulting parent sample contains 144,490 stars with high quality elemental abundances and phase space information, 6,350 of which are located within R GC = 4 kpc of the Galactic centre.

Our first task is to distinguish stars within our parent sample that were formed in situ from those that were likely accreted. For that purpose we display the parent sample on the [Mg/Mn] vs. [Al/Fe] plane in Figure 1 . Das et al. (2020) showed that stars in the APOGEE DR14 sample occupied three major distinct loci on this plane (see also Hawkins et al. 2015) . Most of the stars occupy two large concentrations which, in the DR16 data, are centered at ([Al/Fe],[Mg/Mn]) ∼ (0,0) and (0.3,0.4), corresponding to the low-and high-α disk populations, respectively (e.g., Bovy et al. 2012; Bensby et al. 2014; Nidever et al. 2014; Hayden et al. 2015; Mackereth et al. 2017; Queiroz et al. 2020) . A third group of stars populates a more diffuse "blob" apparent in Figure 1 , which is centered roughly at ([Al/Fe],[Mg/Mn]) = (-0.2,0.5), which Das et al. (2020) associate with an accreted origin. That association was originally proposed by Hawkins et al. (2015) , who identified accreted halo stars as those with high radial velocities, low [α/Fe], and which depart from the relation between radial velocity and Galactic longitude defined by disk stars. The solid line in Figure 1 defines the locus occupied by accreted stars in that chemical plane (top left corner). According to Das et al. (2020) , all stars to the right and below that line are thus deemed to have formed in situ. The latter population is further divided into two sub-classes characterized by low and high α-element abundances, as indicated by the dotted line in Figure 1 . In total our in situ and accreted samples contain 141,514 and 2,976 stars, respectively. The red points on top of the 2D histogram in Figure 1 include only stars within 4 kpc of the Galactic centre. A sizable sub-sample of these so-called bulge stars are located in the accreted locus of chemical space, suggesting that the bulge may host an accreted population. Our sample contains 463 supposedly accreted stars within 4 kpc of the Galactic centre, modulo a contamination by stars formed in situ, which we estimate by different means in Sections 3.1 and 4.2.

Before proceeding with our analysis, it is interesting to examine the association of abundance patterns and stellar origin on the basis of standard Galactic chemical evolution modelling. In this way we hope to ground the empirical definition of accreted vs. in situ chemistry on a theoretical basis. Figure 2 displays the bulge population on the same chemical plane as Figure 1 , but now overlaid by two chemical evolution models starting from the same initial chemical composition, calculated using the flexCE package by Andrews et al. (2017) . The red line shows the "fiducial" model calculated by those authors, with parameters chosen to match the properties of the stellar populations in the solar neighbourhood. We adopt it as a proxy for the behaviour of an in situ population. The blue line represents a model calculated by Hasselquist et al. (2020, in prep.) to match the properties of the Gaia-Enceladus/Sausage system, which we therefore choose to mimic the expected behaviour of accreted populations. The parameters adopted for each model are listed in Table 1 . While the fiducial model adopts a standard exponentially decaying inflow law, that for the accreted population follows a dependence of the te −t type. Both models adopt a Kroupa (2001) IMF with stellar mass ranging from 0.1 to 100 M . In both models the distribution of time delays before the occurrence of SN Ia is an exponential with timescale 1.5 Gyr and a minimum delay time of 150 Myr. Tests were performed to verify that the overall mass scaling is not important for chemical evolution. The deciding factors in fact are the ratio between initial and inflow masses, the outflow mass loading factor, and the star formation efficiency.

The filled circles along the model lines indicate the positions at evolutionary time t = 300 Myr, 1 Gyr, and 5 Gyr. The black cross on both curves indicates the point at which the models reach [Fe/H]=-0.8. The models are far from a perfect match to the data, but provide a reasonable quali- Table 1 . Parameters adopted for chemical evolution models in Figure 2 Parameter

In Situ Accreted Initial gass mass 2 × 10 10 M 3 × 10 9 M Inflow mass scale 3.5 × 10 11 M 6 × 10 10 M Outflow mass loading factor 2.5 6 Star formation efficiency 1.5 × 10 −9 yr −1 1.0 × 10 −10 yr −1 Exponential inflow timescale 6.0 Gyr 2.5 Gyr Figure 2 . APOGEE sample for stars with R GC < 4 kpc compared with chemical evolution models. The red line shows a model for a Milky Way like galaxy and the blue line the model for a dwarf satellite. The models start from the same chemical composition. The filled circles indicate the positions along the evolution after t = 0.3, 1, and 5 Gyr after the beginning of chemical evolution. The black cross in both curves indicates the point at which [Fe/H]=-0.8. Evolution of the massive galaxy is a lot faster, but the star formation rate is much higher, so the two systems are expected to produce similar total stellar mass within the "Accreted" area of the diagram.

tative description of the main trends, which suffices for our purposes.

The first point to be taken from this model comparison is that the "Accreted" region of chemical space in Figure 1 is inhabited by old metal-poor stars from both models, suggesting that both in situ and accreted stars share that area of chemical space. In fact, the accreted stars in halo samples are located in that particular locus of Figure 1 not due to an intrinsic evolutionary property of dwarf galaxies, but rather because their star formation was quenched at the moment when they were disrupted while merging into the Milky Way. As a result, all the stars ever formed by such accreted systems inhabit the unevolved region of Figure 1 . The natural implication of this conclusion is that the "Accreted" region of the chemical plane in Figure 1 is likely to also contain stars that were in fact formed in situ.

Secondly, one can try to exploit the different timescales and star formation rates of the two models to attempt a zeroth order estimate of the expected contribution of in situ and accreted populations to our sample within the "Accreted" region of chemical space. It is clear that chemical enrichment proceeds at a much faster pace in the in situ model, whose gas leaves the "Accreted" region of the chemical plane at t ∼ 250 Myr. Conversely, the gas in the low mass galaxy remains within the "Accreted" locus for about ∼3 Gyr. On the other hand, the star formation rate at the early stages of the evolution of the high mass galaxy, according to the models, is of the order of ∼ 20 M yr −1 , whereas in the case of the dwarf galaxy it is ∼ 0.3 M yr −1 . Therefore, considering the characteristic star formation rates and timescales for the two types of systems, one would expect that the total stellar mass associated with in situ and accreted populations in the "Accreted" region of Figure 1 to be of the same order of magnitude. In Section 6.1.1 we provide a more accurate estimate of the contribution of in situ and accreted populations to the stellar mass budget of the metal-poor bulge, based on their distributions in the E − L z plane.

We conclude that, on the basis of position on the [Mg/Mn] vs. [Al/Fe] plane, it is in principle impossible to completely separate in situ from accreted stars. However, we argue that the clumpy nature of the distribution of stellar populations in Figure 1 speaks in favour of a different origin for at least part of the stellar populations in the "Accreted" region of that plane, and those in the in situ region. In Section 4.2 we show that by imposing selection criteria based on orbital parameters it is possible to keep the contamination of our accreted sample to an acceptable level.

In this section we examine the orbital properties of accreted and in situ stellar populations, as defined in the previous section. We begin this exercise by focusing on the orbital properties of metal-poor populations, with an eye towards establishing whether the chemically defined accreted and in situ groups also differ in their orbital properties. For that purpose, we must first define where we draw the line between metal-poor and metal-rich populations. We begin by examining the metallicity distribution of the sub-populations defined in Figure 1 in order to determine the [Fe/H] threshold defining "metal-poor" stars for the purposes of our exercise. The raw metallicity distribution functions of the accreted, high-and low-α populations are displayed in Figure 3 , together with that for the entire R GC < 4 kpc population. We stress that this raw MDF is not corrected for selection effects, so it is presented here just as a rough guide to the relative metallicities covered by the populations defined in Figure 1 . The MDFs of the high-and low-α populations cover the thin disk, bar, and thick disk structures identified in previous works (e.g., Ness et al. 2013a; Rojas-Arriagada et al. 2017; García Pérez et al. 2018 , Rojas-Arriagada et al. 2020 . The accreted population overlaps substantially in [Fe/H] with the highα population, extending all the way to [Fe/H]∼ −0.8 with a tail towards higher metallicity. As a compromise between obtaining a good representation of the accreted and highα populations, without pushing too far into the metal-rich regime, we adopt an [Fe/H]=-0.8 cutoff. Henceforth, when referring to the "metal-poor" bulge, unless otherwise noticed, we mean stars with [Fe/H]<-0.8. In Section 5.2 we go beyond a simple MDF approach in order to better explore the chemical complexity of bulge stellar populations. Since our goal in this Section is simply to contrast orbital properties of accreted and in situ stars within the same metallicity regime, a straight [Fe/H] cutoff should suffice.

In Figure Of the latter, only 179 belong to the low-α group, as expected given the MDF of the low-α disk (e.g., Hayden et al. 2015) .

Arrows of different colours are placed on the positions of the main in situ and accreted structures identifiable in the E − L z plane, which stand out in the top and bottom panels, respectively. The positioning of the arrows in both panels are exactly the same, to guide the eye. In the top panel, the main in situ structures identifiable are the low-α/thin disk, indicated by the orange arrow, high-α/thick disk stars are indicated by the black arrow. An extended prograde "branch" at about E/10 5 ∼ −1.75 km 2 s −2 , which seems to be associated with the "Splash" population identified by Belokurov et al. (2019) , is indicated by the green arrow. Because of the focus on the metal-poor end of the parent sample, these structures are relatively weak, but in the Appendix the loci of these populations on the IoM place are shown in better detail.

The clumpy distribution of accreted populations in IoM space is clear even under a casual visual inspection of the bottom panel of Figure 4 . By far the most prominent substructure is the large "blob" centered around L z ∼ 0 km s −1 kpc and E/10 5 ∼ −1.5km 2 s −2 (blue arrow) which is associated with the Gaia-Enceladus/Sausage (GE/S) system Belokurov et al. 2018; Helmi et al. 2018; Mackereth et al. 2019) , and is completely absent in the top panel where in situ populations are displayed. Other visually discernible substructures in the bottom panel which are absent in the top panel are the following: (i) another smaller less extended structure to the left of GE/S (E/10 5 ∼ −1.75 km 2 s −2 , and L z /10 3 ∼ −0.6 kpc km s −1 ), which was originally identified by Koppelman et al. (2019) , who named it Thamnos; (ii) A more difuse cloud of prograde stars (L z /10 3 > ∼ 0.6 kpc km s −1 ) at high energy (E/10 5 > ∼ −1.5 km 2 s −2 ), which is associated with the Helmi stream (Helmi et al. 1999) . (iii) a low energy cluster centred at L z /10 3 ∼0.1 km s −1 kpc (red arrow) with energies in the interval −2.6 < ∼ E/10 5 < ∼ −2.0 km 2 s −2 , whose position coincides with that of the low energy (L-E) family of globular clusters identified by Massari et al. (2019) . The bulk of the stars associated with this group are located within ∼4 kpc of the Galactic centre, so we refer to it as the Inner Galaxy Structure (IGS).

One can readily spot significant differences between the two panels, which at face value validate the chemical distinction between accreted and in situ populations. Features associated with well known in situ populations, such as the high-and low-α disks and the Splash (c.f. Figure A2) show much more prominently in the top panel, and nearly vanish in the bottom panel. Conversely, the GE/S system stands out in the bottom panel and is absent in the top panel.

An energy gap is clearly visible in the bottom panel at about E/10 5 ∼ −1.9 km 2 s −2 , which separates stars belonging to the GE/S and the Splash from stars belonging to the IGS. In the top panel, that gap is filled by stars belonging to the thick disk and the Splash, which present a predominantly smooth distribution, much as expected from numerical simulations of Milky Way-like disks undergoing satellite accretion (e.g., Jean-Baptiste et al. 2017 ). This energy gap is not entirely empty in the accreted sample, which is partly due to measurement errors and partly due to the fact that in situ structures, in particular the high-α disk and the Splash make a small, but noticeable, contribution to the data in the bottom panel. This is to be expected on the basis of our discussion in Section 3, where we showed that standard chemical evolution models predict in situ contamination in our accreted sample. Under this interpretation of the data, the distribution of the chemically-defined accreted sample in E − L z space consists of two overlapping populations: one of them consists of various clumps associated with accreted structures, most prominently the GE/S and the IGS, and the other is comprised of a predominantly smooth, metal-poor, residual contamination associated with high-α disk stars, which include both the unperturbed disk and the "Splash".

These results suggest that the IGS is dynamically detached from other metal-poor populations in the halo and disk. This is an important result, so it is vital that the reality of this energy gap is firmed up quantitatively. To do that we measure the ratio between the number of stars with energies within the gap and those below it in both the accreted and in situ samples. We consider stars within the gap to have energies given by −2.00 < E/10 5 < −1.85 km 2 s −2 and those below it have E/10 5 < −2.00 km 2 s −2 . The ratio in the in situ sample is 0.40±0.04 and in the accreted sample 0.27±0.04, configuring a difference that is significant at the 2σ level. We thus conclude that the gap is real, and an indication of the presence of real substructure in the distribution of accreted stars in the low energy locus of IoM space.

Finally, it is worth mentioning that a relatively small number of stars can be seen scattered at high energy levels (E/10 5 > ∼ −1.6 km 2 s −2 ) among the in situ population (top panel), on a wide range of angular momenta. These are predominantly low-α metal-poor stars and are likely accreted contaminants in the in situ sample (c.f. Figure 2 ). . Subsamples used to estimate the contamination of the IGS sample by in situ stars. In both panels, red symbols mark stars located within the locus in orbital space defined by the IGS stars. Blue symbols indicate stars in a reference region of orbital space occupied by high-α stars alone. The ratio between numbers of red and blue stars in the in situ population (top panel) is used to infer the number of in situ contaminants among the IGS stars (red symbols in the bottom panel). We estimate that somewhere between 24 and 40% of stars in our IGS sample are actually highα disk contaminants.

As pointed out in the previous section, examination of Fig Those are metal-poor in situ stars that inhabit the "Accreted" locus of the chemical plane in Figure 1 , as discussed in Section 3.1. We can use those in situ features that are present in the distribution of accreted stars in IoM as a means to estimate the contamination of in situ stars in our IGS sample. For that purpose, we must first establish criteria for IGS membership. We define as IGS members stars meeting the following simple orbital energy and eccentricity criteria:

The energy criterion restricts the sample to the clump that is visually identified in Figure 4 . The eccentricity criterion is aimed at minimising contamination by stars belonging to the inner high-α disk (see Figure A2 ). Next we assess that contamination, using sub-samples of the chemically defined "Accreted" and in situ samples, displayed in Figure 5 , overlaid on the same data originally displayed in Figure 4 . We assume that the energy distribu-tion of stars in the high-α disk is smooth (as predicted by simulations, see, e.g., Jean-Baptiste et al. 2017) and can be estimated from the distribution of in situ stars in the E − L z plane (upper panel of Figure 5 ). We further assume that that energy distribution is the same as that of the contaminants in the accreted sample (bottom panel of Figure 5 ). If those assumptions are correct, one can calculate the ratio in the in situ sample between the numbers of low energy stars in the IGS region, n(low E) InS (red symbols in the top panel of Figure 5 ), and those in a reference high energy locus that is free from accreted stars, n(high E) InS (blue symbols in the top panel of Figure 5 ). That ratio can then be used to convert the number of high energy disk stars in the accreted sample, n(high E) Acc (blue symbols in the bottom panel of Figure 5 ), into the number of contaminating disk stars in the IGS sample, N cont , as follows:

That number can then be compared to the actual number of stars meeting the above defined criteria for the IGS in the accreted region (red symbols in the bottom panel of Figure 5 ) to estimate the contamination fraction.

For that purpose we must select a region of the E − L z space that provides a pure sample of in situ stars belonging predominantly to the high-α disk. We choose such a locus by picking stars with eccentricity < 0.4 and in the energy range −1.9 < E/10 5 < −1.4 km 2 s −2 . Stars in this reference region and in the IGS locus defined above are displayed as blue and red dots in Figure 5 .

By proceeding in this fashion, we obtain the numbers listed in Table 2 , which mean that approximately 22% of our IGS sample actually consist of high-α disk stars. A fundamental systematic uncertainty in this procedure concerns the metallicity cut adopted for the accreted and in situ samples. The calculation above was based on the [Fe/H]=-0.8 definition of metal-poor stars adopted throughout this paper. However, if a slightly more conservative cut is adopted at, say, [Fe/H]=-1.0, the inferred contamination goes up to 40%. This is due to how the numbers in the in situ population, particularly the reference disk stars, depend on metallicity. The angular momentum distribution of disk stars in our sample is strongly dependent on metallicity. At lower metallicity, the L z distribution is shifted towards lower values, leading to a higher [n(lowE)/n(high)] In situ , and thus a higher number of estimated contaminants (N cont ).

We therefore conclude that the contamination of our IGS sample by high-α disk stars ranges between roughly 22 and 40%.

The bulk of the stars belonging to the substructure identified in Section 4.2 are contained within 4 kpc of the Galactic centre. Therefore, it is useful to inspect the orbital behaviour of in situ and accreted stars within that central volume of the Galaxy, with an eye towards determining whether the properties of the accreted population mimic those of their Figure 6 . Distribution of low-α (left panels), high-α (middle panels), and accreted stars (right panels) in action space. Only stars within R GC < 4 kpc are shown.

in situ counterparts, or whether they are distinct and cannot be derived from the latter. Figure 6 displays high-and lowα in situ stars of all metallicities, as well as their accreted counterparts on various IoM spaces. Plotted as a function of angular momentum (L z ) are the vertical action (J z , top panels), radial action (J r , middle panels), and total energy (E, bottom panels). The shift to a focus on the inner 4 kpc causes obvious changes in the distribution of stellar populations. The clearest one is the near complete absence of stars with L z /10 3 > ∼ 1 kpc km s −1 in the inner Galaxy, which is understandable on purely geometric grounds. Interestingly, Splash populations are nearly absent at R GC < 4 kpc. Of the three main in situ structures apparent in Figure 4 , one can only find a mix of the low L z high-α populations and the extension of the low-α towards lower L z . The former also manifests itself by the presence of stars with large J z in the top middle panel. From the low-α disk through the high-α disk to the accreted population, there is a consistent trend towards a larger fraction of the stars having larger vertical and radial actions, as well as lower angular momentum. Accreted populations are dominated by hot, low L z kinematics.

The distribution of stars in the accreted population on all planes differs markedly from that of their in situ counterparts. To gain further confidence on the reality of this critical result, we display in Figure 7 histograms comparing the distributions of action variables of the accreted and in situ inner Galaxy populations. We limit the comparison to the high-α population only, since the differences between accreted and low-α stars can be easily recognised by eye on Figure 6 . Substantial differences can be seen in the distributions of all integrals of motion. This is particulary obvious when comparing the 1D distributions of each integral of motion separately, as shown in Figure 7 . We perform Kolmogorov-Smirnoff (KS) tests comparing the distributions of accreted and in situ populations in those variables, finding that the two samples are not extracted from the same parent population in L z (p-value=1 × 10 −34 ), J z (6 × 10 −15 ) J r (7 × 10 −7 ), and E (5 × 10 −6 ). This result persists, with decreased statistical significance, when only metal-poor stars are considered, with p-value<0.01 in all cases. To isolate stars that are only associated with the IGS, we repeat the KS tests for L z , J z , and J r , by restricting the comparison only to stars with E < −2×10 5 km 2 s −2 . The results are unchanged except for the case of J r , for which the KS test cannot rule out the null hypothesis.

We conclude that the distribution of accreted stars in the inner Galaxy differs from those of major in situ components in action space. This result agrees with the findings by Ness et al. (2013b) and Arentsen et al. (2020b) , who show that metal-poor populations are characterised by hot kinematics (see also Minniti 1996) , even though without detailed chemical compositions they were unable to distinguish in situ from accreted populations. In the most detailed investigation of the kinematics of bulge stellar populations to date, Ness et al. (2013b) found that the metal-poor component of the bulge ([Fe/H] < ∼ -1) does not partake in its cylindrical rotation, being characterised by slightly lower rotation and higher velocity dispersion (see also Zasowski et al. 2016) .

We point out that the differences between accreted and high-α in situ populations are reduced when only metalpoor ([Fe/H]<-0.8), low-energy stars (E/10 5 < −2km 2 s −2 ) are considered. In that case, the only variable for which a residual statistically significant difference is found is L z (p-value=0.004). This is an interesting result. It suggests that, if indeed the IGS is the remnant of an accretion event, it has had the time to mix well with its co-spatial in situ population, which suggests a very early accretion event. We discuss this and other possible interpretations of the data in Section 6.1 .

In summary, the chemically defined accreted populations in the bulge are not dynamically associated with their in situ counterparts. They are dynamically hotter than the in situ population typically possessing low angular momen- Figure 7 . Comparison between the distributions of high-α in situ and accreted inner Galaxy populations in action variables. Statistically significant differences can be seen in L z and J z , whereby accreted stars are dynamically hotter, being distributed towards smaller angular momentum (L z ) and larger vertical action (J z ). The energy distributions differ substantially due to the presence of GE/S stars with E> −2 × 10 5 km 2 s −2 . Removal of that contamination does not change the results. tum and being distributed towards large values of vertical and radial action. We take these results as suggestive of an accretion origin. In the next section we discuss the chemical properties of the IGS.

In the previous section we investigated the orbital properties of the IGS, concluding that they are distinct from those of co-spatial in situ populations, which supports an accretion Figure 9 . Comparison between stars from the IGS (red dots) with their bulge counterparts (gray dots). The IGS has an abundance pattern that is characteristic of those of dwarf satellites of the Milky Way, yet they seem to form a single sequence with the remainder of the bulge populations, even though they all are associated with different Galactic components.

origin. In this section we examine the chemical compositions of the stars from the IGS alongside those of their inner Galaxy counterparts. Our main objective is to check whether detailed chemistry is consistent or not with an accretion scenario. In Figure 8 we display again the accreted stars shown in the bottom panel of Figure 4 , this time marking stars deemed to be members of the IGS and GE/S systems as red and blue dots, respectively. The IGS sample is defined in Section 4.3. We consider members of the GE/S system to be those stars with −0.4 < L z /10 3 < 0.3 kpc km s −1 and −1.75 < E/10 5 < −1.3 km 2 s −2 . According to these criteria, our samples contain 259 candidate members of the IGS and 1,026 associated with the GE/S system. These definitions are henceforth adopted in our discussion of the chemical properties of the IGS.

In Figure 9 we contrast the stars from the IGS with their in situ R GC < 4 kpc counterparts in chemical space. We limit our comparison to a few elements that sample a range of nucleosynthetic pathways, with Mg being mostly contributed by SN II, C and N by massive and AGB stars, and Ni by SN Ia (e.g., Chiappini et al. 2003; Nomoto et al. 2013) . The latter element is additionally interesting due to it being typically depressed in dwarf satellites of the Milky Way (e.g., Shetrone et al. 2003) . For more comprehensive studies of the detailed abundance patterns of bulge stars, we refer the reader to the excellent studies by Schultheis et al. (2017) and Zasowski et al. (2019) and the review by Barbuy et al. (2018) .

Two main results emerge from this comparison: (i) the stars belonging to the IGS have an abundance pattern that is typical of dwarf galaxies, characterised by low abundances of elements such as Al (by construction), C, N, and Ni (e.g., Hayes et al. 2018; Mackereth et al. 2019; Helmi 2020; Das et al. 2020) ); (ii) the abundance pattern of the IGS apparently behaves as a natural extension of the bulge field population towards [Fe/H] < −1, a result that has also been noted by Schultheis et al. (2017) and Zasowski et al. (2019) . To be clear, this extension of the trends towards metal-poor populations concerns only the locus occupied by stars in the various chemical planes, not the relative counts. The latter are examined in Section 5.2. Despite this apparent continuity, it is worth noting that slight deviations of this pattern are apparent in elements such as Mg, Ni, and Si (not shown). Similar behaviour is seen in the data presented by Queiroz et al. (2020) 4 , who examine the spatial variation of the stellar distribution on the α-Fe plane across the Galactic disk and bulge. Naively, one would consider this seeming continuity between the IGS and its in situ counterparts surprising, given the alleged accreted nature of the IGS. In Section 5.1 we show that this seeming continuity between the chemical compositions of the IGS and its more metal-rich bulge counterparts is actually expected on the basis of numerical simulations.

It is also worth pointing out that IGS stars with [Fe/H]>-1.0 show a variety of different behaviours. Some of them follow the same trend as the high-α disk in Mg, Ni, and C+N, whereas others display a behaviour more similar to those of dwarf satellites, with low [Mg/Fe], [Ni/Fe], and [(C+N)/Fe]. By construction, all of them display low [Al/Fe]. We hypothesise that these metal-rich IGS stars are a mixture of thick disk contaminants (whose low [Al/Fe] ratios are difficult to understand), thin disk contaminants, and genuine IGS members.

In summary, the results presented in this Section and Sections 4.2 and 4.4 show that a metal-poor stellar population located in the inner Galaxy shares chemical and IoM properties with well known halo accreted systems. Figure 9 suggests on the other hand that the abundance pattern of this same structure is consistent with sharing a common chemical evolution path with its bulge counterparts potentially associated with other co-spatial Galactic components. Yet the analysis in Sections 4.2 and 4.4 suggests that the IGS is dynamically detached from the rest of the bulge (Figures 6  and 7 ). This is a puzzling result. In the next subsections we look into the theoretical expectations for the chemical evolution of Milky Way-like galaxies that underwent comparable accretion events (Section 5.1) and inspect the data for the presence of a chemical connection between the accreted and in situ populations in the inner Galaxy.

The discussion above naturally prompts us to enquire into the theoretical predictions for the chemical evolution of Milky Way-like Galaxies with accretion histories similar to that seemingly implied by our results. For that purpose we resort to the predictions from the suite of cosmological numerical simulations by the Evolution and Assembly of GaLaxies and their Environments (EAGLE) project (Schaye et al. 2015; Crain et al. 2015) . The question guiding this analysis is the following: what are the chemical properties of Milky Way like galaxies undergoing an important early accretion event? More specifically, do we expect the chemical properties of the stars belonging to the accreted system to differ substantially from those of the in situ population?

The EAGLE simulations provide excellent material for this enquiry. The large volumes included in the various sets of simulations make possible robust statistical analyses of the evolutionary outcomes for a broad gamut of final galaxy properties, and at a wide enough range of resolutions to address the needs of different types of investigations. To maximise resolution we choose to study the 25 3 cMpc 3 L025N752-Recal simulation, whose volume is large enough to yield a substantial number of galaxies with Milky Waylike masses and accretion histories similar to that suggested by our results. The chemical evolution prescriptions of the simulations track only elements relevant to the thermal balance of the interstellar medium, yet this is sufficient for our purposes, as predictions for crucial elements Fe, C, N, and Mg are included. For more details on the simulations we refer the reader to Schaye et al. (2015) and Crain et al. (2015) .

The analysis was conducted as follows. We selected galaxies with Milky Way-like masses within the interval 7 × 10 11 < M 200 /M < 3 × 10 12 and looked for those that underwent an accretion of a satellite with stellar mass M /M ≥ 10 8 at z ≥ 2 (see Section 6.2 for an estimate of the mass of the IGS progenitor). We found a total of 15 accreted galaxies with such a profile (34% of the host galaxy sample within that halo mass range). The particles for the accreted galaxies were selected when it was at its peak stellar mass. For three of those galaxies we display chemical composition data in Figure 10 . The data are displayed in the form of contour maps indicating the distribution of stars on abundance planes, with in situ and accreted populations shown in gray and red, respectively. In situ populations in this context are defined as stellar particles that were bound to a halo in the main galaxy branch at the snapshot prior to star formation. The abundance planes shown are Mg-Fe, C-Fe, and N-Fe, which track the contribution of SN II, SN Ia, and AGB stars, sampling a range of chemical enrichment timescales.

The answer to the question posed above immediately strikes the eye upon inspection of Figure 10 . For all three galaxies, and all three elemental abundances, the abundance ratios of the accreted populations are in good agreement with those of the in situ populations at same metallicity. This is the case for all the 15 galaxies meeting our selection criteria. It is also the case for all the other elemental abundances tracked by EAGLE. This result is also fully consistent with the analytical model presented in Figure 2 , which shows that the early chemical evolution of a Milky Way-like galaxy and its satellite are similar. We conclude therefore that the seeming continuity between the abundances of the accreted and in situ populations in Figure 9 is in fact not surprising, but actually an expectation from theory, at least for the elements whose chemical evolutions are modelled by EAGLE. This is a very important result, the implications of which configure the proverbial double-edged sword. On one hand the theoretical predictions for the locus of accreted populations in chemical space cannot rule out the accreted nature of the IGS. On the other, it cannot ascertain it either, as the apparent continuity between the chemistry of the IGS and its in situ populations could alternatively be a by-product of pure in situ chemical evolution.

To find our way out of this conundrum we invoke relative number counts in chemical space. In trying to understand the apparent chemical continuity between the accreted and in situ populations in Figure 9 we ask the following question: is the distribution of stellar data on abundance planes consistent with the accreted population being an additional structure, or do the number counts suggest that it is just the tail end of the distribution of the in situ populations? In the next section we address this question in a quantitative fashion.

In this Section we tackle the problem of understanding the apparent chemical continuity between the IGS and its more metal-rich inner Galaxy counterparts ( Figure 9 ) from a purely observational perspective. We center our approach around the following question: can the stellar density distribution in chemical space be described by a discrete number of detached substructures? Ness et al. (2013a) showed that the bulge MDF can be well modelled by a combination of five components, which they ascribe, in order of increasing metallicity, to the inner halo, the metal-poor thick disk, the bulk of the thick disk itself, the bar, and the thin disk. Some of these peaks can be identified in our raw MDF in Figure 3 MDF is quite difficult due to differences in the two metallicity scales and in the relative strengths of the various peaks (recall that our MDF is not corrected for selection effects), and the uncertainties associated with multiple gaussian fitting in 1D space. These difficulties are particularly important towards low metallicity, where the statistics is not as robust and the results of abundance analysis may be afflicted by important systematics.

The exceptional quality of our detailed chemical compositions comes in handy in this juncture. By considering the distribution of our sample in higher dimensional space we hope to gain additional leverage in distinguishing the different components making up the stellar population mix of the inner Galaxy. The stellar distribution on some chemical planes, such as [Mg/Mn]-[Al/Fe] and [Al/Fe]-[Fe/H], show a clumpy distribution which suggests that the data may be amenable to such an approach. In other words, one should be able to model the stellar density distribution in chemical space as a combination of overlapping structures, which are not necessarily connected to one another. To test this hypothesis we check whether the stellar density distribution in a few chemical planes can be reproduced by a Gaussian Mixture Modelling (GMM) approach. Obviously a Gaussian distribution is only a zero-th order approximation to the stellar density distributions in chemical space. A massive stellar system such as, for instance, the Galactic disk, underwent star formation long enough that the run of the number of stars with metallicity and abundance ratios is strongly determined by the star formation history. Conversely, in a less massive system whose star formation was quenched during an early accretion to a massive halo, the distribution of stellar metallicities should be narrower and thus the GMM approach may work better.

We are particulary concerned with modelling the metalpoor end of the metallicity distribution to check whether the IGS is detached from its more metal-rich peers in the inner Galaxy. The key measure of chemical detachment of the IGS from its more metal-rich counterparts in the inner Galaxy is provided by the ability of a GMM approach to account for the density distribution of the IGS, so that no important residual connection between that and the other populations is left after model subtraction (although some residual should be expected given the contamination by in situ populations discussed in Section 3).

The sample for analysis is displayed in the left panels of Figure 11 . It consists of all stars within R GC < 4 kpc, but removing outliers with [Al/Fe]>+0.5, as their abundances are affected by the N-rich-star phenomenon (Schiavon et al. 2017) . We also removed relatively metal-poor stars with low Figure 11 ). On the other planes, accreted stars are located in the low metallicity end of the distribution ([Fe/H] < ∼ -1). We fit a separate GMM to the distribution of these stars on each of the three chemical planes. The ideal number of components was decided by calculating the Bayesian Information Criterion (BIC) for each fit. We found that the best match is found when either 4 or 5 components are adopted. We allow the mean and covariances of each component to be entirely free with no informative or uninformative priors. The GMMs for the 5 component case are overlaid on the data in the middle panels of Figure 11 . The standard deviation of the residuals around the fits are 0.23, 0.19, and 0.17 for the top, middle, and bottom panels, respectively.

The residuals after subtraction of the best fit models and the data are displayed in the right panels. The residuals are relatively small, and slightly larger towards higher metallicity, which reflects the much larger stellar density in that regime. When normalized by the Poisson fluctuations per bin, the residuals are homogeneously distributed across the planes, and are everywhere smaller than ∼ 2σ of the Poisson uncertainty in each bin. These numbers suggest that the distribution of the stellar density in the locus of accreted stars is well accounted for by the GMM approach, showing that they can be matched by a simple combination of detached, yet overlapping structures. Interestingly, the GMM fit to the stellar distribution in the Mg-Fe plane highlights the fact that the lowest metallicity component has a slightly lower [Mg/Fe] than the next one in order of increasing metallicity. This result supports the notion that the most metal-poor component is dominated by an accreted stellar population, as one would expect [Mg/Fe] to be constant in an in situ population within the metallicity interval involved. Also quite importantly, even though the GMM match was performed independently to the data on the three chemistry planes, the associations of stars to a given group were largely consistent across the three models, which attests for the reality of these groups.

The five-component match identifies peaks located at [Fe/H] ∼ -1.3, -0.8, -0.5, and -0.2 and +0.3, which are very roughly consistent with the five MDF peaks identified by Ness et al. (2013a) , after allowing for a ∼0.2-0.4 dex difference between the zero points of the metallicity scales of the APOGEE and ARGOS surveys (see, for instance, the positions of the peaks in the APOGEE bulge MDF determined by García Pérez et al. 2018, and Rojas-Arriagada et al. 2020, in prep.) . The bulge MDF has been analysed exhaustively in recent work based, e.g., on the ARGOS (Ness et al. 2013a) , Gaia-ESO (Rojas-Arriagada et al. 2017) , and APOGEE (Schultheis et al. 2017; García Pérez et al. 2018) surveys. The number of peaks in the bulge MDF, which in turn constrains the number of physically distinct components cohabiting the bulge, have been subject to debate in the recent literature (see discussion by Fragkoudi et al. 2018 ). We do not aim at a detailed analysis of the topic in this paper. However, we wish to highlight the fact that the consideration of additional abundances helps untangling the two most metal-poor populations in our sample. The task of disentangling these populations in the MDF of Figure 3 is far from straightforward, as the distributions of these components overlap strongly in one-dimensional metallicity space. However, consideration of additional elemental abundances lifts the degeneracy between these components, bringing them to sharp relief on elemental abundance planes such as those in the middle and bottom rows in Figure 11 .

The analysis presented in Section 5.2 leads to a few important conclusions: (i) the distribution of the numbers of stars over various chemical planes can be modelled by the combination of overlapping yet detached sub-structures. This is particularly true in the low metallicity end, where residuals are the lowest, adding further support to the notion that the IGS has an accreted origin; (ii) the fit in the Mg-Fe suggests that the [Mg/Fe] distribution of the lowest metallicity component is suggestive of an accretion origin; and (iii) the GMM method consistently picks out the same substructures in various chemical planes, involving the abundances of Fe, Mg, Al, and Mn. The metallicities of these substructures are a reasonably good match to those of well known peaks in the MDF of the Galactic bulge (Ness et al. 2013a ). This latter result validates the application of the GMM method to describe the distribution of stellar data in chemical space.

We start this section by briefly summarising our results before proceeding with a discussion of their implications. In Sections 4 and 5 we discussed evidence for the presence of a remnant of an accretion event in the inner Galaxy, on the basis of the existence of substructure in orbital and chemical space. While the whole body of evidence for accretion is strong, when considered in isolation, our interpretation of either dynamical or chemical evidence may reasonably be called into question.

The evidence from the orbital side hinges on two main results: (i) the statistically significant differences in orbital properties between the accreted stars and their more metalrich in situ counterparts in the inner Galaxy (Figures 6 and  7) ; (ii) the existence of an energy gap separating the IGS from other metal-poor halo substructure. The presence of this energy gap seems to be immune to selection effects, since it is absent in chemically defined metal-poor in situ populations (Figure 4) , whose spatial distribution, uncorrected for selection effects, is the same as that of their accreted counterparts.

It is important to note, however, that when R GC < 4 kpc stellar populations are considered in isolation, there is a smooth dependence of orbital properties on metallicity, which could be accommodated by a fully in situ scenario. On the other hand such smooth behaviour would also be expected in the accreted case, given that the orbital and chemical properties of substructures overlap strongly due to a combination of observational error and cosmic variance.

On the chemistry side, evidence for accretion relies on (i) the similarity between the chemical composition of the IGS and those of satellites of the Milky Way and other accreted systems (Figure 13 ), and (ii) the presence of clumpiness of the stellar distribution in various chemical spaces (Figures 1 and 11) . Regarding (i), we showed that the chemistry of unevolved in situ populations are also similar to that of their accreted counterparts (Figure 2) . Moreover, there is an apparent continuity between the chemical compositions of the IGS and its bulge counterparts (Figure 9 ). However we also showed that similarity between the chemical compositions of accreted and in situ populations is an actual expectation from cosmological numerical simulations (Figure 10) . Regarding (ii) we showed that the data in chemical space can be well matched by the combination of detached substructures ( Figure 11 ). However, one cannot rule out that such clumpiness is the result of a bursty history of star formation and chemical enrichment.

Despite the above caveats, when considered in aggregate, chemical and orbital information suggest the presence of a previously unidentified structure located in inner Galaxy, the IGS. We hypothesise that this structure is the remnant of a galaxy accreted to the Milky Way in its very early life. In the next Subsections we discuss the properties of the hypothetical progenitor of the IGS. In Section 6.1 we discuss the origin of the metal-poor populations in the Galactic bulge, in Section 6.2 we discuss the properties of the putative progenitor of the IGS; and in Section 6.3 we present a speculative scenario for the early evolution of the Milky Way halo.

In this Section we estimate the contribution of accreted and in situ stars to the metal-poor star content of the Galactic bulge, and contrast that estimate with predictions from fully in situ chemical evolution models. Figure 1 . The 2D histogram shows metal-poor stars at R GC > 4 kpc, and the red (green) dots represent accreted (in situ) stars in the inner Galaxy. A clear bimodal distribution in [Al/Fe] is visible in the data for both the inner and outer Galaxy. Right panel: Same stars displayed on the energy-angular momentum plane. The data for the outer Galaxy are dominated by the GE/S and Splash systems. In the inner Galaxy, the accreted sample is dominated by the IGS, but about 1/10 of the sample is consistent with an association with the GE/S system (E > −1.85 × 10 5 km 2 s −2 ). The in situ sample is dominated by thick disk stars, although a small fraction seems to belong to the Splash. Note that the distributions of accreted and in situ stars in IoM space are different at a statistically significant level (see Section 4.2).

We begin by assessing the relative contribution of accreted and in situ stars to the metal-poor stellar mass budget of the inner Galaxy. The data shown in Figure The right panel of Figure 12 shows the distribution of accreted and in situ stars on the E − L z plane. The accreted star sample is dominated by low energy stars, mostly associated to the IGS. About 10% of the accreted stars are located at high energies and are likely stars associated to the GE/S system (large gray blob) which are currently crossing the inner Galaxy. Notably, very few in situ stars (green points) have such high energies. They predominantly occupy the low energy region of the plane, with a number of stars positioned along the Splash (to the right of GE/S, see green arrow in Figure 4 ).

Other studies have spotted duplexity in the properties of metal-poor bulge samples, which may be attributed to the presence of a mix of accreted and in situ populations. Pietrukowicz et al. (2015) identified two major sequences in the period-amplitude plane for a large sample of fundamental-mode RR Lyr variable stars. The two populations differ slightly in metallicity, with the more metalpoor population outnumbering their metal-rich counterparts by approximately a factor of 4. Those authors suggest that this duality is indicative of merging activity in the early history of the Milky Way bulge. Kunder et al. (2020) obtained radial velocities for a few thousand RR Lyrae stars from the Pietrukowicz et al. (2015) sample and determined the presence of two kinematically distinct populations. One population follows bar-like orbits and the other one is more centrally concentrated, being characterised by hotter kinematics. They suggest that the latter may result from an accretion event that preceded the formation of the bar.

On the theory side, cosmological numerical simulations predict relatively intense accretion activity in the first few Gyr of massive galaxies' lives. They have reached a degree of sophistication that enables a quantitative confrontation with our data. Early numerical simulations by Kobayashi & Nakasato (2011) propose that the bulge was fully formed by the merger of galaxies with M ∼ 5 × 10 9 M at z > ∼ 3. On the other hand, the more recent EAGLE simulations provide more detailed information that can be compared with our data. For instance, in the simulated galaxies whose abundances are displayed in Figure 10 , the fraction of accreted stars within 4 kpc of the centre for [Fe/H] < −1 varies from 0.10 to 0.55. Along the same lines, the analysis of FIRE simulations by El-Badry et al. (2018) finds that same fraction within 1 kpc of simulated Milky Way-like galaxy centres to be around 0.5. More recently, Fragkoudi et al. (2020) analysed high resolution data for Milky Way-like galaxies from the Auriga suite of cosmological numerical simulations. They found that the fraction of the stellar mass with [Fe/H]<-1 within 4 kpc of the centre that has an accretion origin to vary between 0.13 and 0.8. Moreover, for the two simulated galaxies matching the properties of the Milky Way best, that fraction was smaller than 0.4.

Considering all the uncertainties on the observational and theoretical fronts, and taking into account the relative arbitrariness of the spatial and [Fe/H] limits adopted in such comparisons, our results broadly confirm the predictions by numerical simulations from various groups. Obviously, better statistics on both the observational and numerical simulation fronts will be required for such data to provide stringent constraints on models.

It is encouraging that the main tentative result presented in this paper does not fall into a theoretical void, that it is in consonance with expectations from cosmological numerical simulations. Recourse to theoretical predictions to validate an observational result is nonetheless methodologically, and indeed philosophically, unsound. Thus, before discussing the properties of the putative progenitor of the IGS, it behooves us to mention theoretical predictions by various groups which match data for bulge stellar populations without resort to satellite accretion.

The history of modelling the chemo-kinematic properties of bulge stellar populations is decades long, so that paying tribute to all the excellent work done in the past is beyond the scope of this paper. For that we refer the reader to Section 4 of the review by Barbuy et al. (2018) . We restrict our discussion to recent work by Matteucci et al. (2019 Matteucci et al. ( , 2020 , because to our knowledge it is the only instance of a detailed comparison betwen model predictions and multiple chemical compositions for a large sample based on APOGEE data. Their analytical model for the history of star formation and chemical enrichment of the bulge stellar populations provides a good match to the MDF (Matteucci et al. 2019 ) and the abundances of O, Mg, Si, Ca, Al, K, Mn, Cr, and Ni. While such models resort to ad hoc assumptions about gas infall for different Galactic components, and some non-negligible amount of yield adjustment, they match the observed data very well. The implication is that the fraction of the stellar mass within the bulge that is due to accretion is zero at all metallicities. In that scenario, the entirety of the bulge stellar populations, including those deemed accreted by our criteria described in Figure 1 , were actually formed in situ.

In light of our discussion in Section 3, where we showed that in situ and accreted populations share a similar path in chemical space (Figure 2) , the results by Matteucci et al. (2019 Matteucci et al. ( , 2020 are not surprising. As pointed out in Section 6, orbital information such as that discussed in Section 4 is indeed required for that distinction to be made. Moreover, recall that according to our criteria there are 463 accreted stars out of a total sample of 6,350 stars within 4 kpc of the Galactic centre. Assuming a contamination at the 22-40% level of the accreted sample by in situ thick disk stars, the total accreted stellar mass amounts to only ∼4-5% of our sample 5 . Therefore, as far as the chemical evolution of the bulge is concerned, Matteucci's work offers, at face value and without consideration of the Milky Way's insertion into the cosmological context, a viable in situ formation scenario that can possibly explain the chemical compositions of vast majority of the stellar populations therein. Notwithstanding this remarkable success, accounting for the considerable morphological and kinematic differences between the overlapping components within a few kpc of the Galactic centre under a single unified chemodynamical evolutionary picture remains a challenge.

While the IGS contributes a very small fraction of the stellar mass of the Galactic bulge, it is a major contributor to the stellar mass budget of the stellar halo. Therefore, it is important to give further consideration to the reality of its discrimination as a separate system, which is the topic of the next sub-section.

Recent studies based on N-body simulations (e.g., Jean-Baptiste et al. 2017; Koppelman et al. 2020 ) have shown that substructure in kinematic and orbital planes cannot be ascribed to individual accretion events in a straightforward manner. In particular, Jean-Baptiste et al. (2017) simulated the evolution of orbital properties of accreted and in situ populations in the event of accretion by single and multiple satellites with varying merger mass ratio. The main conclusions from that work have consequences for the interpretation of our results.

Firstly, Jean-Baptiste et al. (2017) show that distributions on the E − L z plane that resemble qualitatively those presented in Figures 4 and 8 can result from the accretion of one single satellite with a 1:10 mass ratio. It can be seen from their Figure 2 that several Gyr after the occurrence of the accretion event, satellite remnants are distributed across a large range of energies and angular momenta, spanning the loci occupied by both the IGS and GE/S in Figure 8 . Most of the substructure in orbital space in their simulations is associated to the accreted populations, including low energy clumps resembling that characterising the IGS in Figure 4 . By the same token, the distribution of the in situ population is a lot smoother and strikingly similar to the distribution in the upper panel of Figure 4 (minus the Splash). Their simulations also show that accretion moves stars formed in situ into high energy and low angular momentum orbits, where accreted stars are usually found. In fact, in situ stars contribute to some of the substructure on the E − L z plane, although not at the same level as the accreted populations, which dominate the substrucutre, particularly at low energy. Along the same lines, Koppelman et al. (2020) analysed a set of N-body simulations that match the observed properties of the GE/S debris today. According to those simulations, the core of the accreted satellite spirals in and is dissolved, with some of the stars acquiring slightly prograde orbits.

On the basis of such simulations one might reasonably argue in favour of a scenario in which the GE/S and the IGS are both part of the same progenitor. In that situation, one could conceive that the core of GE/S, being more massive and denser than the outskirts, could have sunk into the inner Galaxy under the effect of dynamical friction, leaving its less bound stellar tenants behind to be observed in higher energy orbits today, where GE/S stars are found (c.f. Figure 4) .

It is not easy to reconcile such a scenario with the relative distribution of IGS and GE/S data in chemical space in Figure 13 . Stars belonging to the IGS on average present higher mean [Mg/Fe] and lower mean [Fe/H] than their GE/S counterparts. As pointed out in Section 6.2, the IGS lacks a substantial "knee-shin" component to its distribution in the Mg-Fe plane, unlike the GE/S. Note that these differences exist despite the fact that the chemical criteria defining both samples are identical. The IGS is clearly less chemically evolved than the GE/S. Therefore, if the IGS was indeed the core of the GE/S, that hypothetical progenitor system would be characterised by a positive metallicity gradient in the early universe. Galaxies in the nearby and intermediate redshift universe are predominantly characterised by negative metallicity gradients (e.g., Stott et al. 2014; Goddard et al. 2017) . That is also the behaviour of satellites of the Milky Way such as the Magellanic Clouds (e.g., Cioni 2009 ) and the Sagittarius Dwarf (e.g., Hayes et al. 2020) . This is understood as being due to the fact that star formation in galaxies takes place in an inside-out manner, since their central regions evolve more quickly than their outskirts, due to availability of larger amounts of gas for star formation. On the other hand, higher redshift galaxies with positive metallicity gradient have been identified, an occurrence that may be caused by strong winds associated with intense star formation in the centres of dwarf galaxies at z ∼ 2 (e.g., Wang et al. 2019) . It is not clear however how prevalent that phenomenon is.

Secondly, it is worth highlighting additional results presented by Jean-Baptiste et al. (2017) where simulated Milky Way-like galaxies undergo 1:10 mergers with more than one satellite. Those lead to very thorough mixing of the stars from different satellites in orbital space, particularly in the low energy region of the E − L z plane, where the IGS was identified. This suggests that our accreted sample should include stars from more than one satellite in the low energy regime. Under that light our result above is surprising. Even if the IGS and the GE/S were indeed separate systems, we would expect that some of our low energy accreted stars show an abundance pattern consistent with an association with the GE/S. On the contrary, it is only the high energy stars (E/10 5 > ∼ −2 km 2 s −2 ) in our bulge accreted sample (c.f. Figure 12 ) that possess abundances that are consistent with an association with the GE/S. We conclude that it is unlikely, though not impossible, that the IGS is the core of the Gaia Enceladus/Sausage system. That would require that the progenitor of the Gaia Enceladus/Sausage system have a slightly positive metallicity gradient and have a core that is less chemically evolved than its outskirts, which is uncommon.

In the previous Sections we showed that the inner ∼4 kpc of the Milky Way hosts a stellar population with the following properties: (i) it presents a chemical composition that resembles those of accreted systems (e.g., Hayes et al. 2018; Mackereth et al. 2019; Helmi 2020) ; (ii) it is detected as substructure in IoM space; (iii) on average it presents slower rotation, higher vertical action, and slightly higher radial action than metal-rich bulge populations; (iv) it is chemically detached from its metal-rich counterparts. In aggregate these properties suggest the presence of an accreted stellar population within 4 kpc of the Galactic centre. In this Section we use the observed properties of the remnant population to estimate the mass and history of the putative progenitor system. This is achieved through comparison with data for a better known accreted system (GE/S) and the predictions of the EAGLE simulations.

In Figure 13 we compare the chemical properties of the IGS with those of the GE/S system Belokurov et al. 2018; Helmi et al. 2018; Mackereth et al. 2019 ), a massive satellite (M ∼ 10 8.5 M Mackereth & Bovy 2020) , that was accreted to the Milky Way about 8-10 Gyr ago. For details on the selection of GE/S stars, see Figure 8 and related discussion. We ignore, for the sake of simplicity, any selection effects that may bias the relative distribution of chemical properties of the observed stellar populations associated with the two galaxies. In other words, we assume that the samples displayed in Figure 8 represent the bulk of the populations of the two systems. Consequently the conclusions drawn from comparison of the IGS and GE/S should be considered qualitative, and subject to revision should that assumption be demonstrated to be incorrect. With that caveat in mind, Figure 13 displays stars from the two systems on various abundance planes. The IGS is represented by red symbols, the GE/S as blue symbols, and the parent sample, which is dominated by stars from the high-and low-α disks, is plotted as gray-scale 2D histograms. The data show that stars belonging to the IGS have slightly higher average Mg, Ni, C, and N abundances than those of GE/S at [Fe/H] ∼ −1. This result suggests that the progenitor galaxy of the IGS was likely more massive than GE/S. We recall that Mackereth et al. (2019) used predictions from the EAGLE numerical simulations to estimate the stellar mass of GE/S. In particular, they showed that the mean [Fe/H] and [Mg/Fe] of the stellar particles associated with dwarf simulated galaxies were good indicators of stellar mass (see Figure 11 of Mackereth et al. 2019) .

We resort to that predicted relation between mass and chemistry in order to estimate the mass of the progenitor of the IGS. We used the scipy curve_fit optimisation routine to fit the data for satellites of Milky Way-like galaxies in the L0025N075-RECAL simulation (see Mackereth et al. 2019 

Where the zero-th order term has been adjusted to match the mass inferred by is likely to have been an important contributor to the stellar mass budget of the Milky Way halo.

The distribution of the data for the IGS in the Mg-Fe plane may be used to constrain the history of star formation of the progenitor, and thus the time of accretion. One cannot help but notice the fact that the data for the IGS and GE/S differ in one fundamental aspect: our IGS sample has almost no stars with [Fe/H] > ∼ −1, and some of them are likely thick disk contaminants (see Section 5). That is the regime where the GE/S displays a strong decline in [Mg/Fe] for increasing [Fe/H], associated with the onset of enrichment by SN Ia, which takes place on a timescale of ∼ 10 8−9 year (e.g., Maoz et al. 2012; Nomoto et al. 2013 ). Assuming we are not being deceived by an unforeseen selection effect, the absence of this more metal-rich population in the data for the IGS may indicate an early quenching of star formation, likely associated with the accretion to the Milky Way at a very early time.

This result confirms theoretical expectations from the E-MOSAICS simulations by Pfeffer et al. (2020) and the EAGLE predictions for the chemical compositions of early accreted satellites ( Figure 10 ). Based on an analysis of the expected orbits of globular clusters connected with accreted satellites in the simulations, those authors concluded that the orbits of surviving accreted clusters tend to have lower energies if their host satellites were more massive and/or were accreted at higher redshifts. Furthermore, this observed truncation at [Fe/H] ∼ −1 is in agreement with results from a previous study by our group (Horta et al. 2020) , where the chemical compositions of Galactic GCs in APOGEE were studied. In that paper we showed that the more metal-poor half of the GCs associated with the L-E subgroup, whose position in IoM space coincides with that of the IGS, have chemical compositions that are consistent with an accreted origin. It is possible that the more metal-rich GCs have an in situ origin, associated with the thick disk (see Figure 3 from Horta et al. 2020 ).

As pointed out in Section 4, the locus of the IGS in IoM space matches closely that of the L-E family of globular clusters identified by Massari et al. (2019) . Those authors showed that the L-E clusters followed a tight age-metallicity relation, suggesting that they might be associated with a massive satellite. Indeed, both the orbital properties and time of accretion of the progenitor of the IGS agree with the predictions by Kruijssen et al. (2019b Kruijssen et al. ( , 2020 for the Kraken accretion event. According to those authors, the Kraken was a satellite of the Milky Way that was accreted within the first few Gyr of the Galaxy's existence (see also Forbes 2020) . The prediction is based on a detailed comparison of the ages, metallicities, and orbital properties of the Milky Way globular cluster system with predictions of the E-MOSAICS suite of cosmological numerical simulations (Pfeffer et al. 2018; Kruijssen et al. 2019a) . Although the specific identities of the globular clusters associated to the Kraken or the L-E group are not unique, their median metallicities are roughly similar to that of the IGS ([Fe/H]∼ −1.5 and −1.24, respectively). The orbital properties of those globular clus- Figure 13 . Comparison between the IGS (red symbols) and Gaia-Enceladus/Sausage (GE/S, blue symbols) in chemical space. The data for the parent sample are also plotted (gray), to mark the positions of the low and high-α disks in these diagrams. ters, with apocentric distances between 3 and 5 kpc, are also consistent with an association with the IGS.

Finally, Kruijssen et al. (2020) estimate the stellar mass of the Kraken to be ∼ 2 × 10 8 M , or approximately 70% of the mass of GE/S. Our estimate from equation (1) is that the progenitor of the IGS is about twice as massive as the Kraken, and approximately twice the mass of GE/S. Considering the uncertainties associated with the various methods, we do not deem the absolute differences worrying. Further details on the physical properties of the progenitor of the IGS, as well as a definitive association with the Kraken will have to await a more careful examination of both upcoming new data and predictions from improved numerical simulations.

With the above caveats in mind, we conclude by summarising our speculative views on the most likely story implied by the data presented in this paper. While the vast majority of the stars in the Milky Way were likely formed in situ (e.g., Bland-Hawthorn & Gerhard 2016), it is generally agreed that, under a Λ-CDM framework, galaxies such as the Milky Way experienced intense accretion activity in their early lives (e.g., Naab & Ostriker 2017) .

In Section 6.1.1 we found that somewhere between 1/4 and 1/3 of the metal-poor stars 7 in the inner Galaxy were accreted, and the rest were formed in situ. Our inability to perfectly distinguish one group from the other leaves a key question in the air: is there an evolutionary connection between the two or were they formed and evolved in entire isolation before co-habiting the inner few kpc of the Galaxy? A possible sequence of events, is one whereby the most metal-poor stars in the inner Galaxy were formed in situ, followed by the early accretion of the progenitor of the IGS which, under the effect of dynamical friction, took a quick dive towards the heart of the young Milky Way. Assuming the IGS brought with it fresh gas 8 , it is likely that dissipation brought it to the Galactic centre, triggering a burst of star formation which would influence the later course of the chemical evolution in the inner Galaxy, explaining the chemical continuity apparent in Figure 9 9 .

It will be interesting to see whether such speculative scenarios will survive scrutiny by more detailed observations and sophisticated modeling and, if so, how it will jibe with our understanding of the formation, evolution, and interaction between other structures co-habiting the inner Galaxy.

We report evidence to the existence of an accreted structure located in the inner Galaxy (the IGS). To reveal its existence we resorted to a search for clustering in multidimensional parameter space, based on the best available data from APOGEE and Gaia, and astroNN distances. This structure consists of a metal-poor population with detailed chemistry resembling that of low mass satellites of the Milky Way with low abundances of Al, C, N, and Ni. It is chemically and dynamically detached from the more metal-rich populations with which it shares its location in the heart of the Galaxy, and from other accreted populations across the Milky Way halo. We suggest that this population is the remnant of a satellite that was accreted early in the history of the Milky Way. The chemistry of this population is characterised by a plateau of relatively high [Mg/Fe] and the absence of a "knee" towards higher [Fe/H] indicating an early quenching of star formation activity, which is consistent with a very early accretion event. Comparison with mean chemical composition predictions from numerical simulations, we infer a stellar mass of M ∼ 5 × 10 8 M , or approximately twice that of the Gaia-Enceladus/Sausage system. If indeed its existence is confirmed by further investigation, one would 7 Which themselves amount to only ∼5% of the stellar mass budget of the Galactic bulge (see Ness et al. 2013b) 8 This assumes ram pressure stripping was not important, which may be reasonable. Ram pressure stripping is proportional to v 2 and since the IGS is located at the Galactic centre, the putative merger must have occurred at a low velocity. 9 Along similar lines, Whitten et al. (2019) proposed that the steep age gradient of stellar populations in the inner halo (R GC < 14 kpc), inferred from BHB star mean color gradients, can be possibly associated with occurrence of a few dissipative mergers in the early Milky Way history. be led to conclude that its progenitor was a major building block of the Galactic halo, which has until now remained elusive due to its location behind a curtain of extinction and foreground crowding. The orbital and chemical properties of the IGS, the implied time of accretion, and the mass of its progenitor are in good qualitative agreement with those predicted by Pfeffer et al. (2020) and Kruijssen et al. (2020) . We summarise below the main implications and questions raised by this finding:

• According to our favoured interpretation of the data, the IGS is a remnant of a satellite galaxy accreted to the Milky Way a long time ago. That such an accretion event happened in the Milky Way's very early history is an expectation from cosmological numerical simulations that predict an intense accretion activity in massive galaxies at z ∼ 2 − 3 (e.g., El-Badry et al. 2018; Pfeffer et al. 2020; Fragkoudi et al. 2020) . The Milky Way itself may have been the subject of unusually strong accretion activity at those early times, which may explain observed properties of its disk and halo populations (e.g., Mackereth et al. 2018; Schiavon et al. 2020 );

• We present evidence for an accretion event that took place before the main structures of what we understand today as the Galactic bulge were ever formed. Thus, the IGS should be viewed as a major building block of the Galactic halo, rather than a very minor component of the Milky Way bulge, to whose stellar mass budget it contributes no more than a few percent. The vast majority of the stars in the bulge today were most likely formed in situ (e.g., Kunder et al. 2012; Di Matteo 2016; Debattista et al. 2017; Fragkoudi et al. 2018) .

• The abundance ratios of stars associated with the IGS are well aligned with the sequence occupied by more metalrich populations located within the inner ∼4 kpc of the Galactic centre. Naively, that would be a surprising result, as one may expect that accreted and in situ populations might have undergone distinct early histories of star formation and chemical enrichment. We show that numerical simulations actually predict such seeming continuity between the loci occuppied by accreted and in situ populations in all abundances tracked by the EAGLE simulations.

• As a logical path for further enquiry, one would like to explore at least three possible avenues: (i) verify whether models can match the distribution of stellar density in chemical space, which may be a key discriminator between accretion and in situ formation scenarios; (ii) contrast accreted and in situ populations in other abundance planes including elements that sample additional nucleosynthetic pathways, such as those resulting from s-and r-process production. (iii) contrast the age distributions of stars in the chemodynamically defined accreted and in situ populations. All these lines of investigation will benefit from additional data made available in forthcoming observational efforts.

• The combined mass of the progenitor of the IGS and the GE/S adds up to 8 × 10 8 M . Considering mass estimates for other accreted systems (e.g., Koppelman et al. 2019; Forbes 2020; Kruijssen et al. 2020 ) and dissolved globular clusters (Schiavon et al. 2017 , Horta et al., 2020 , the total mass in accreted satellites may be well in excess of ∼ 1.6 × 10 9 M . Although all these mass estimates are subject to considerable zero point uncertainties, this result is, at face value, in mild tension with independent estimates of the total stellar mass of the Galactic halo (e.g., .

• There seems to be an interesting connection between the finding reported in this paper and the presence of a population of N-rich stars in the inner Galaxy (Schiavon et al. 2017) , which is generally interpreted as resulting from the destruction of an early population of globular clusters (e.g., Martell & Grebel 2010; Martell et al. 2016; Koch et al. 2019; Savino & Posti 2019; Hanke et al. 2020; Fernández-Trincado et al. 2020) . The frequency of this population has been shown to be larger by an order of magnitude within the inner few kpc of the halo than in its outer regions (Horta et al. 2020, in prep.) . Moreover, Kisku et al. (2020, in prep.) show that roughly 1/4 of the N-rich population within ∼ 4 kpc of the Galactic centre has elemental abundances consistent with an accretion origin, and orbital properties that are similar to those of the IGS. In future work we intend to examine how these finding jibe together visà vis our current understanding of the accretion history of the Milky Way and models for globular cluster formation and destruction in a cosmological context. This paper presents an exploration of the chemical complexity of the metal-poor bulge on the basis of a statistically significant sample. The window into the bulge metal-poor bulge populations will be further widened in the next data release of APOGEE 2. Moreover, the Pristine Inner Galaxy Survey (PIGS) is performing AAOmega+2dF spectroscopy with the 3.9 m Anglo-Australian Telescope of ∼8,000 targets selected on the basis of CaHK photometry collected with MegaCam on the Canada-France-Hawaii Telescope, (Arentsen et al. 2020a ) amassed the largest sample to date of very metal-poor stars (∼1,300 stars with [Fe/H]< −2) in the inner Galaxy. Detailed abundance analysis of such a sample will certainly shed light onto the nature of the early populations of the bulge. Further into the near future, the MOONS 10 surveys of the Milky Way (González et al. 2020, Messenger, in press.) will likely increase the sample of bulge metal-poor stars with detailed chemical compositions by an order of magnitude. Figure A1 . Definition of high-and low-α disks on the Mg-Fe plane. We added a ±0.04 dex "cushion" in [Mg/Fe] around the dividing line in order to minimize inter-contamination between the two samples.

where the lines defining each sub-sample are drawn on top of the data. Metal-poor stars are removed in order to minimize contamination by halo stars.

High-and low-α disk stars are displayed in action and energy space in Figure A2 . The arrows indicate the approximate loci of stars belonging to the each sub-sample. The orange arrow points to the location of low-α disk stars (e.g., Bovy et al. 2012; Bensby et al. 2014; Nidever et al. 2014; Hayden et al. 2015; Mackereth et al. 2017; Queiroz et al. 2020) , which are characterised by high angular momentum, small energy scatter, and low radial (J r ) and vertical (J z ) actions. There are, however, outliers towards higher E, J r , and Jz, which are due to contaminants from the GE/S system and possibly other accreted systems.

The high-α disk is dominated by two major structures. For constant L z , these two high-α subgroups occupy significantly different and well separated values of J r and energy. The most important for the purposes of this paper is strongly rotational, overlapping substantially with the lowα disk at high L z , but also extending towards L z ∼ 0 and presenting larger scatter towards higher E, J r , and J z (magenta arrow). The other component presents higher energy at L z /10 3 < ∼ 1 kpc km s −1 and merges with the other components at higher L z , being associated with the Splash population (see Belokurov et al. 2019 ) (green arrow). This paper has been typeset from a T E X/L A T E X file prepared by the author. Figure A2 . Distribution of high-and low-α disk populations in action space. The orange arrow indicates the rough position of the low-α/thin disk populations, the magenta arrow that of the standard high-α/thick disk populations, and the green arrow that of stars associated with the Splash.

Multi-Object Spectroscopy in the Next Decade: Big Questions, Large Surveys, and Wide Fields

The Messenger

SciPy: Open source scientific tools for Python

NumPy: A guide to NumPy

The Galactic Bulge

IAU Symposium

Astronomical Society of the Pacific Conference Series

We thank the many people around the world whose hard work fighting the ongoing COVID-19 pandemic has made it possible for us to remain safe and heathy for the past several months. (Hunter 2007) , Galpy (Bovy 2015; Mackereth & Bovy 2018), TOP-CAT (Taylor 2005) .

In this appendix we briefly present the distribution of disk stars in action space, to inform the discussion of accreted and in situ populations in Section 4. We start by distinguishing high-and low-α disk stars. For that purpose we display our parent sample on the Mg-Fe plane in Figure A1 ,