# The Pristine survey

## XXIV. The Galactic underdogs: Dynamic tales of a Milky Way metal-poor population

Isaure González Rivera de La Vernhe<fnref n="1" />,*, Vanessa Hill<fnref n="1" />, Georges Kordopatis<fnref n="1" />, Felipe Gran<fnref n="1" />, Emma Fernández-Alvar<fnref n="2" />,<fnref n="3" />,<fnref n="1" />, Anke Ardern-Arentsen<fnref n="4" />, Guillaume F. Thomas<fnref n="3" />,<fnref n="2" />, Federico Sestito<fnref n="5" />, Camila Navarrete<fnref n="1" />, Nicolas F. Martin<fnref n="6" />,<fnref n="7" />, Else Starkenburg<fnref n="8" />, Akshara Viswanathan<fnref n="8" />, Giuseppina Battaglia<fnref n="3" />,<fnref n="2" />, Kim A. Venn<fnref n="5" />, and Sara Vitali<fnref n="9" />,***

Received 25 April 2024 / Accepted 24 June 2024

## ABSTRACT

**Context.** Metal-poor stars hold key information on the early Milky Way. Through the identification and characterisation of substruc- tures, one can understand internal mechanisms (including merger and accretion events), which are indispensable to reconstruct the formation history of the Galaxy.

**Aims.** To allow an investigation of a population of very metal-poor stars ([Fe/H] < -1.7) with disc-like orbits (planar and prograde), high angular momenta ($L_z/J_{tot} > 0.5$) and rotational velocities ($V_\phi > 180 \text{ km.s}^{-1}$) proposed in the literature, we used a sample of ~3 M giant stars with Gaia DR3 BP/RP information and Pristine-Gaia metallicities down to -4.0 dex that we aimed to decontaminate. To achieve this, we constructed a sample as free as possible from spurious photometric estimates, an issue commonly encountered for high $V_\phi$ metal-poor stars.

**Methods.** We created a statistically robust sample of ~36 000 Pristine-Gaia very metal-poor ([Fe/H] < -1.7) giant stars, using APOGEE and LAMOST data (adding GALAH and GSP-spec for verification) to estimate and remove contamination. We investi- gated the spatial and kinematic properties of the decontaminated sample, making use of $V_\phi$ as well as the action space, which are both powerful tools to disentangle stellar populations.

**Results.** The global distribution of very metal-poor stars in our sample shows the typical kinematics, orbital properties, and spatial distributions of a halo; however, as in previous works, we found a pronounced asymmetry in the $L_z$ and $V_\phi$ distributions, in favour of prograde stars. We showed that this excess is predominantly due to prograde-planar stars (10% of the very metal-poor popula- tion), which can be detected down to [Fe/H] = -2.9 at a $2\sigma$ confidence level. This prograde-planar population contains stars with $V_\phi > 180 \text{ km.s}^{-1}$ and $Z_{max} < 1.5 \text{ kpc}$. While the overall orbital configurations ($Z_{max} - R_{max}$ or action space distributions) of our sample match that of a halo, the highly prograde and planar subset (2% of the very metal-poor population) also bears characteristics classically associated with a thick disc: (i) a spatial distribution compatible with a short-scaled thick disc, (ii) a $Z_{max} - R_{max}$ distribution similar to the one expected from the thick disc prediction of the Gaia Universe Model Snapshot, and (iii) a challenge to erase its signature assuming a stationary or prograde halo with $\bar{V}_\phi \sim 30-40 \text{ km.s}^{-1}$. Altogether, these results seem to rule out that these highly prograde and planar stars are part of a thin disc population and, instead, support a contribution from a metal-weak thick disc. Higher resolution spectra are needed to fully disentangle the origin(s) of the population.

**Key words.** catalogs – surveys – Galaxy: disk – Galaxy: formation – Galaxy: kinematics and dynamics

## 1. Introduction

Our galaxy, the Milky Way (MW), is a bridge between near and far-field cosmology. Within the scope of ΛCDM, a hierarchical formation scenario of galaxies can be constrained and compared to high-redshift observations (Freeman & Bland-Hawthorn 2002; Frebel & Norris 2015). In addition, there is an opportu-

nity to unveil the early history of the MW through the study of its metal-poor ([Fe/H] < -1.0) stars. The atmospheres of the most metal-poor stars, reflective of limited chemical enrichment within their birth gas clouds, is the reason why they are com- monly associated with the oldest populations of the Universe (Beers & Christlieb 2005; Frebel & Norris 2013, 2015).

Until recently, it was widely accepted that the most metal-poor (i.e. oldest) stars would be primarily found on pressure-supported orbits (Eggen et al. 1962; Tumlinson 2010)

<footnote n="1">¹ Université Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Laboratoire Lagrange, 06304 Nice, France</footnote>

<footnote n="2">² Departamento de Astrofísica, Universidad de La Laguna, E-38206 La Laguna, Tenerife, Spain</footnote>

<footnote n="3">³ Instituto de Astrofísica de Canarias, E-38200 La Laguna, Tenerife, Spain</footnote>

<footnote n="4">⁴ Institute of Astronomy, University of Cambridge, Madingley Road, Cambridge CB3 0HA, UK</footnote>

<footnote n="5">⁵ Department of Physics and Astronomy, University of Victoria, PO Box 3055, STN CSC, Victoria BC V8W 3P6, Canada</footnote>

<footnote n="6">⁶ Université de Strasbourg, CNRS, Observatoire astronomique de Strasbourg, UMR 7550, F-67000 Strasbourg, France</footnote>

<footnote n="7">⁷ Max-Planck-Institut für Astronomie, Königstuhl 17, D-69117 Heidelberg, Germany</footnote>

<footnote n="8">⁸ Kapteyn Astronomical Institute, University of Groningen, Lantleven 12, 9747 AD Groningen, The Netherlands</footnote>

<footnote n="9">⁹ Instituto de Estudios Astrofísicos, Facultad de Ingeniería y Ciencias, Universidad Diego Portales, Av. Ejército Libertador 441, Santiago, Chile</footnote>

<footnote>* Corresponding author; isaure.gonzalez-rivera@oca.eu</footnote>

<footnote>** Millenium Nucleus ERIS</footnote>
– either within the Galactic halo or the bulge (White & Springel 2000; Starkenburg et al. 2017a; El-Badry et al. 2018). This idea follows the expectations of a hierarchical assembly process (Searle & Zinn 1978; Schuster & Nissen 1989; Ibata et al. 1994; Chiba & Yoshii 1998; Zolotov et al. 2009; Belokurov & Kravtsov 2022; Xiang & Rix 2022), and also highlights the contribution of the halo to a metal-poor inner Galaxy (Arentsen et al. 2020; Lucey et al. 2021; Ardern-Arentsen et al. 2024). The search for metal-poor stars in our galaxy goes back to the 1950s (see Beers & Christlieb 2005 for a complete review). The pioneering work of Schwarzschild & Schwarzschild (1950) and Roman (1950) paved the way for the discovery of the first very metal-poor (VMP<fnref n="1" />) stars ([Fe/H] < -2.0, Bond 1970, 1980). However, most works from that epoch emphasise the rare nature of this population, making their investigation a rather tedious task in terms of sample selection, [Fe/H] validation, and abundance determination. During the late 1990s to early 2000s, the development of digital photometric and spectroscopic instruments marked the beginning of the era of large ground-based surveys dedicated to the search of metal-poor stars. These include the HK survey (Beers et al. 1985, 1992) and the Hamburg/ESO survey (HES, Wisotzki et al. 1996, 2000), which are responsible for the large samples of stars being uncovered with even lower metallicities, namely extremely metal-poor (EMP) stars ([Fe/H] < -3.0). In particular, HES is at the origin of the class of ultra metal-poor (UMP) stars ([Fe/H] < -4.0) with the discovery of HE 0107-5240 (Christlieb et al. 2004).

The mid-2010s represent a turning point in galactic archaeology. With the advent of the space-based Gaia mission (Perryman et al. 2001; Gaia Collaboration 2016a,b), one could get easy access to the 3D positions of millions of celestial objects within ~20 kpc and with unprecedented accuracy, leading to the serendipitous detection of metal-poor, VMP, and EMP stars. Gaia's third data release (Gaia Collaboration 2023a) exploits the medium resolution spectra (R ~ 11 500) obtained by the Radial Velocity Spectrometer (RVS, Cropper et al. 2018; Katz et al. 2023) to provide one of the largest samples of stellar parameters and elemental abundances available (Recio-Blanco et al. 2016, 2023; Gaia Collaboration 2023b). This combination of astrometric and chemical parameters has offered new insights into the identification and validation of the most metal-poor stars of the Galaxy, and eventually on their chemodynamical characterisation (Matsuno et al. 2024; Viswanathan et al. 2024a, see also Venn et al. 2020; Kielty et al. 2021). In parallel, this database can be complemented with recent ground-based spectroscopic surveys, whether high resolution, such as APOGEE DR17 (R ~ 22 500, Majewski et al. 2017; Abdurro'uf et al. 2022), GALAH DR3 (R ~ 28 000, De Silva et al. 2015; Buder et al. 2021), and Gaia-ESO (R ~ 19 800–21 500 with GIRAFFE and R ~ 47 000 with UVES for MW targets, Randich et al. 2022), or low resolution, such as RAVE DR6 (R ~ 7500, Steinmetz et al. 2020) and LAMOST (R ~ 1800, Cui et al. 2012; Chen et al. 2022). As part of the Sloan Digital Sky Survey (SDSS, Abazajian et al. 2003), we can also cite SEGUE (R ~ 1800, Lee et al. 2008), BOSS, and eBOSS (R ~ 1560–2270 in the blue channel, Dawson et al. 2013; Comparat et al. 2016). Consequential progress in the identification of metal-poor stars is being made thanks to the development of photometric surveys such as SkyMapper (Keller et al. 2007), S-PLUS (Mendes de Oliveira et al. 2019), and more specifically

with the Pristine survey, a photometric survey at the CFHT using a narrow-band filter centred on the CaHK doublet (λ ~ 3900 Å).

Interestingly, several authors are now challenging the formation history of the Galactic disc (Sestito et al. 2019; Di Matteo et al. 2020; Mardini et al. 2022; Bellazzini et al. 2024). It has been established for a long time that the disc is a combination of two main components (Gilmore & Reid 1983): the thin disc (scale height ~ 300 pc), with circular orbits and a metallicity distribution peaking close to solar values (Hayden et al. 2014, 2015; Buder et al. 2019; Mackereth et al. 2019), and the thick disc, vertically more extended (~1000 pc) with slightly lower metallicities down to ~1.5 dex (Norris et al. 1985; Morrison et al. 1990; Chiba & Beers 2000; Ruchti et al. 2011; Kordopatis et al. 2011, 2013). Very recently, several studies identified a significant number of metal-poor stars with kinematics typical of the thin and thick discs (Beers et al. 2014; Sestito et al. 2019, 2024; Di Matteo et al. 2020; Cordoni et al. 2021; Zhang et al. 2023; Dovgal et al. 2024; Nepal et al. 2024; Gallart et al. 2024), in particular through the eyes of Pristine (Sestito et al. 2020; Venn et al. 2020; Fernández-Alvar et al. 2021, 2024; Viswanathan et al. 2024b). This population, as well as its prograde nature, is supported by recent cosmological simulations (Santistevan et al. 2021; Sestito et al. 2021; Sotillo-Ramos et al. 2023; Carollo et al. 2023), but its full characteristics are not yet well established.

However, one would expect very ancient stars rotating in a disc to have heated up towards halo-like kinematics over time, and not to be found on disc-like orbits. A first explanation for such a divergence between observations and theory is that the sample of Pristine low-metallicity candidates is contaminated by a fraction of the overwhelmingly larger population of metal-rich stars ([Fe/H] > -0.5). In fact, the metallicity distribution function (MDF) of the MW in different Galactic locations shows a dominant metal-rich population related to the thin disc (Hayden et al. 2015; Kordopatis et al. 2015). Inevitably, the contamination of metal-poor samples by metal-rich counterparts is higher at rotational velocity ($V_φ$) values typical of the thin disc. Indeed, metal-poor samples with prograde and cold kinematics (small velocity dispersion), that is, disc-like, are considerably more polluted by metal-rich stars compared to metal-poor samples located at higher latitude $b$ and, or with kinematically hotter orbits. Such an effect, adding to the fact that extinction biases on Pristine photometric metallicities are larger along the disc area, represents a serious concern, as misclassifying only a very small fraction of the metal-rich sample could compromise the detection of a kinematically cold metal-poor population. Leaving aside doubts related to the quality of photometric metallicities, the origin of this disc-like metal-poor population is still highly debated. Some authors relate it to a primordial disc, with divergent opinions on its in situ or ex situ origin (Sestito et al. 2021; Mardini et al. 2022; Bellazzini et al. 2024), some to a prograde halo (Di Matteo et al. 2020; Belokurov & Kravtsov 2022; Zhang et al. 2023), while others involve the Galactic bar in the migration of metal-poor stars from the bulge or the halo towards the disc (Yuan et al. 2024; Li et al. 2024a; Dillamore et al. 2023).

In this paper, we aim at characterising a population of VMP stars confined on prograde-planar orbits, specifically with high angular momenta and high $V_φ$, similarly to MW disc stars. Using the Pristine-Gaia synthetic photometric metallicity catalogue (PGS), which is part of the first Pristine data release (Martin et al. 2024), we conscientiously selected a sample of about 3M giant stars with available radial velocities (RVs) from Gaia (Katz et al. 2023), in the range -4.0 < [Fe/H] < 0. The main emphasis of this

<footnote n="1">¹ In the following sections, we refer to VMP stars as stars with [Fe/H] < -1.7.</footnote>
work is to obtain a sample that is statistically significant on the metal-poor end to confidently assess the existence of a VMP disc population, while being as free as possible from metal-rich outliers at high $V_φ$. This step is the most important to be able to discuss the potential origins of this peculiar population.

This paper is organised as follows. Section 2 presents the relevant catalogues needed to construct a sample with astrometric and chemodynamical information. Section 3 describes the quality cuts and the filtering method applied to the sample, and Sect. 4 details the validation of the latter. Section 5 is dedicated to the decontamination of the sample using spectroscopic surveys. Finally, in Sects. 6 and 7 we discuss and summarise our results.

## 2. The data: PGS metallicities and Gaia RVS orbits

The Pristine survey (Starkenburg et al. 2017b) is a narrow-band photometric survey, observed with the Canada-France-Hawaii Telescope (CFHT), and focused on the metallicity-sensitive Ca H & K doublet lines (located at 3933.7 and 3968.5 Å), which are of great interest to infer photometric metallicities down to [Fe/H] ~ -4.0 dex. This combination of narrow-band and broad-band photometry has led to the successful discovery of large samples of EMP, VMP and UMP stars in recent years (Starkenburg et al. 2018; Aguado et al. 2019; Bonifacio et al. 2019; Venn et al. 2020; Lardo et al. 2021; Kielty et al. 2021; Lucchesi et al. 2022; Caffau et al. 2023; Viswanathan et al. 2024a).

For this work, we used the publicly available PGS metallicity catalogue, an all-sky metallicity catalogue based on Gaia DR3 BP/RP spectro-photometry published by Martin et al. (2024) as part of Pristine DR1. Relevant columns for this paper are listed in Table 1. In this catalogue, 52 million Gaia DR3 objects with BP/RP coefficients<fnref n="2" /> have their de-reddened synthetic CaHK magnitudes, CaHK_0, computed after integrating their BP/RP spectra under the transmission curve of the same CaHK filter as the one used in the Pristine survey. Methods are described in Martin et al. (2024) following Gaia Collaboration (2023c) for the integration of photometric filters, and in Starkenburg et al. (2017b) for the metallicity estimation.

Compared to Pristine CFHT observations (6500 deg<fnref n="2" /> mapped since 2015), PGS allows for a complete sky coverage: the catalogue is a great opportunity to gather a statistically significant sample of stars with a net prograde motion, close to the Galactic plane, while displaying very low metal contents. Such a population was highlighted in Fernández-Alvar et al. (2021) from CFHT Pristine observations. In order to confirm these results, especially at low metallicities, we dedicated a large part of this work to the preliminary verification of the catalogue parameterisation.

To be able to investigate both the kinematical (velocities) and dynamical (actions) characteristics of this population, we complemented our sample with orbital parameters derived using Gaia RVs and the geometric distances of Bailer-Jones et al. (2021) as our distance estimator. Galactocentric cartesian positions X, Y, Z, Galactocentric cylindrical radius R, and Galactocentric cylindrical velocities $V_R$, $V_φ$, $V_Z$ were computed for all stars with available RVs. Orbital parameters were obtained following Kordopatis et al. (2023), who used the Stäckel fudge method for their computations (Binney 2012a; Bovy & Rix

2013; Mackereth & Bovy 2018), along with the Galpy<fnref n="3" /> code (Bovy 2015). The adopted MW potential is the McMillan (2017); it is axisymmetric, with a modelling of cold gas discs, a thin and a thick disc, a bar, a bulge, and a halo. The potential was adapted to recent values of Solar positions, that is, $(R, Z)_⊙ = (8.249, 0.0208)$ kpc (GRAVITY Collaboration 2020; Bennett & Bovy 2019), and Solar velocities, that is, $(V_R, V_φ, V_Z)_⊙ = (-9.5, 250.7, 8.56)$ km.s⁻¹ (GRAVITY Collaboration 2020; Reid & Brunthaler 2020). Further details on the computation are provided in Sect. 3.1 of Kordopatis et al. (2023).

The following sections describe the different quality cuts applied. The number of stars removed from each individual cut is summarised in Table A.1.

## 3. The sample: Data selection

### 3.1. Gaia DR3 quality cuts

After cross-matching PGS with the Gaia DR3 source catalogue, a filtering based on astrometry has been applied:

- $\bar{\omega} > 0$ mas to avoid negative parallaxes<fnref n="4" />;

- relative parallax error, that is, $\delta\bar{\omega}/\bar{\omega} < 0.20$ to avoid large errors propagated onto velocities<fnref n="5" />;

- sqrt(phot_g_nobs)/phot_g_meanaim flux over error < 0.015 to minimise the presence of variable stars;

- Renormalised Unit Weight Error (RUWE) < 1.4 to avoid astrometric binaries, following Lindegren et al. (2021a);

- we removed stars with source_id in gaiadr3.vari_summary, Gaia DR3's variable sources catalogue (complementary to Pristine's Pvar, see Table 1).

### 3.2. Pristine DR1 quality cuts

We adopted the following cuts relative to the PGS catalogue (see also Sect. 7.3 of Martin et al. (2024) and Table 1):

- mcfrac_CaHKsyn > 0.9 to remove stars with a large fraction of their Monte Carlo samplings falling beyond the edges of the resulting metallicity grid (Martin et al. 2024 advise mcfrac_CaHKsyn > 0.8, we adopted a stricter cut since 80% of the catalogue's mcfrac_CaHKsyn lies above 0.9);

- FeH_CaHKsyn_50th > -4.0 to avoid stars at the lower edge of the metallicity grid adopted for the training sample;

- Pvar < 0.3 to avoid photometrically variable stars<fnref n="6" />;

- d_CaHK < 0.05 mag to remove noisy photometric measurements;

- Schlegel et al. (1998) extinction value $E(B-V) < 0.5$ mag<fnref n="7" /> following Sect. 7.5 of Martin et al. (2024).

<footnote n="3">³ http://github.com/jobovy/galpy</footnote>

<footnote n="4">⁴ Quality cuts relative to parallax do not take into account the correction from the parallax zero point of -17 µas (Lindegren et al. 2021b).</footnote>

<footnote n="5">⁵ When applying a stricter cut on the parallax (e.g. $\delta\bar{\omega}/\bar{\omega} < 0.10$), we restrict the accessible volume, therefore our sample selection sharply decreases (by 24%), which hampers our ability to probe the population of disc-like metal-poor stars we focus on. However, the conclusions of this work remain unchanged.</footnote>

<footnote n="6">⁶ As discussed in Sect. 5 of Martin et al. (2024), the variability of stars is complex on the bright end of the G band, when obtained from Gaia data. It is necessary to develop a model that is dependent on the noise budget of the source to determine which Gaia stars are intrinsically variable.</footnote>

<footnote n="7">⁷ Applying a stricter cut on $E(B-V)$ (e.g. <0.2) results in depleting the final metal-poor sample statistically, but in a homogeneous way, preserving the distributions in kinematics discussed in the following sections.</footnote>

<footnote n="2">² In the Gaia archive, BP/RP spectra are stored as coefficients with their corresponding covariance matrix (projections on a set of basis functions), which can be reconstructed as spectral energy distributions by a linear combination.</footnote>
<table><caption>**Table 1.** Description of key columns from the published *Pristine-Gaia* synthetic photometric metallicity catalogue.</caption><thead><tr><th>Column</th><th>Description</th><th>Unit</th></tr></thead><tbody><tr><td>G_0</td><td>De-reddened G magnitude</td><td>mag</td></tr><tr><td>BP_0</td><td>De-reddened BP magnitude</td><td>mag</td></tr><tr><td>RP_0</td><td>De-reddened RP magnitude</td><td>mag</td></tr><tr><td>CaHK_0</td><td>Synthetic CaHK magnitude</td><td>mag</td></tr><tr><td>d_CaHK</td><td>Uncertainty on synthetic CaHK magnitude</td><td>mag</td></tr><tr><td>FeH_CaHKsyn</td><td>Photometric metallicity obtained directly from the Pristine model</td><td>dex</td></tr><tr><td>FeH_CaHKsyn_50th</td><td>Photometric metallicity obtained after 100 Monte-Carlo draws of the Pristine model</td><td>dex</td></tr><tr><td>mcfrac_CaHKsyn</td><td>Fraction of Monte-Carlo draws falling into the [Fe/H] grid</td><td>–</td></tr><tr><td>Pvar</td><td>Probability for a source to be variable</td><td>–</td></tr><tr><td>E(B-V)</td><td>Schlegel et al. (1998) extinction value</td><td>mag</td></tr></tbody></table>

*Pristine-Gaia* synthetic provides two metallicity estimates: FeH_CaHKsyn_50th, the median photometric estimate obtained after 100 Monte-Carlo draws of the giant Pristine model<fnref n="8" />, and FeH_CaHKsyn, the direct estimate (see Table 1). We found that, irrespective of the location in the PGS colour-magnitude diagram (CMD, in this context, refers to a colour-absolute magnitude diagram), there is an overall good agreement between median and direct metallicity estimates, with an absolute difference consistently below 0.04 dex, except for a few sources with (BP - RP)₀ < 0.8 and M<sub>G₀</sub><fnref n="9" /> ≥4 mag, for which the difference is greater than 0.1 dex (0.83% of the total sample). For the remainder of this work, we define FeH_CaHKsyn_50th as our generic photometric metallicity, that we dub [Fe/H]PGS, to be compared with spectroscopic metallicities.

### 3.3. Stellar population quality cuts and isochrone filtering method

Once our data have been curated based on Gaia DR3 astrometry (Sect. 3.1) and PGS criteria (Sect. 3.2), additional filters were applied:

- we removed stars identified in recent censuses of globular (Vasiliev & Baumgardt 2021) and open clusters (Hunt & Reffert 2023);

- we selected giant stars using the separation curve in the CMD of Fernández-Alvar et al. (2021) (blue curve in the right panel of Fig. 1). As advised in Starkenburg et al. (2017b) and Martin et al. (2024), the Pristine survey was designed to target mostly metal-poor giant stars of the MW halo and is consequently less precise for main sequence stars (smaller training sample).

An additional filter was created to remove stars with underestimated PGS metallicities, based on their location in the CMD. Indeed, it is expected that genuine metal-poor giants (intrinsically old) occupy the bluest part of the red giant branch (RGB). For that reason, one can use theoretical isochrones in order to define areas in which genuine metal-poor stars (in that context, stars with estimated [Fe/H]PGS < -1.5) should not reside in. However, on the one hand, metal-poor stars can have a wide range of ages; on the other hand, PGS metallicities can have rel-

atively large uncertainties ($\delta_{[Fe/H]}$ ~ 0.3 dex). So, the definition of this area needs to be quite generous. That is why we adopted a PARSEC isochrone<fnref n="10" /> (Bressan et al. 2012) of [Fe/H] = -0.5 and $\tau$ = 5 Gyr to define the physical boundaries into which a PGS star with [Fe/H] < -1.5 can reside (see Fig. 1).

We used an argument based on age and metallicity to justify the physicality of the location of our sample in the CMD. If a star within the metal-poor subsample appears redder than the boundary defined by the isochrone, then either:

- its metallicity is underestimated, and is therefore higher than the metallicity we are targeting for the metal-poor subsample;

- its isochrone age is older than that of the Universe ($\tau$ = 13.5 Gyr);

- the distance or the de-reddening procedure is inaccurate, causing erroneous metallicity and, or orbits.

With this method, we filtered out 15% (9811 sources) of the PGS metal-poor subsample. Eventually, the PGS sample comprises a total of 2880 338 giant stars, including 54 009 metal-poor stars according to our definition. In particular:

- stars with [Fe/H]PGS < -1.5 were filtered according to Sects. 3.1, 3.2 and 3.3;

- stars with [Fe/H]PGS > -1.5 were filtered similarly without the isochrone filtering;

- we note that this sample does contain horizontal branch (HB) stars; however, they only make up for 11% of the total size.

The spatial distribution of our PGS sample is displayed in Fig. 2. It is essentially distributed within the extended Solar neighbourhood, but covers quite well the X-Y space (left panel), and also between 2 < R < 14 kpc and -6 < Z < 6 kpc (middle panel). As seen from the R-Z distribution, most stars lie close to the plane, although due to Pristine extinction cuts, there are very few stars below |Z| = 2 kpc at large R. As expected, we can clearly see the thin and thick disc stand out in the $V_\phi$ - [Fe/H] space (right panel), but we do see a tail of distribution extending down to the lowest metallicities, as discussed in Fernández-Alvar et al. (2021). Additionally, we note the presence of a ‘pencil-beam’ structure in X-Y that could be caused by the Gaia scanning law (see Sect. 5.2 of Gaia Collaboration 2016b), in particular because of our sample selection, which only comprises stars with Gaia RVS RVs (RVS selection effects can also be a cause for the structure). We discuss these asymmetries in Sect. 6.1.

<footnote n="8">⁸ The Pristine model focusses on giant stars, thus the generic output for the photometric metallicities is the one computed with the giant model (FeH_CaHKsyn_50th). However, we note that we chose to assign to dwarf stars the photometric estimate obtained with the dwarf model (FeH_CaHKsyn_dw_50th).</footnote>

<footnote n="9">⁹ We used the geometric distances of Bailer-Jones et al. (2021) to compute M<sub>G₀</sub>.</footnote>

<footnote n="10">¹⁰ http://stev.oapd.inaf.it/cgi-bin/cmd</footnote>
<figure><img src="image_1.png" /><figcaption>**Fig. 1.** Colour – magnitude diagram (CMD) of the *Pristine-Gaia* synthetic (PGS) metallicity catalogue. The sample, after applying the quality cuts described in Sects. 3.1 and 3.2 is colour-coded in grey. *Pristine-Gaia* synthetic stars with [Fe/H]$_{PGS}$ < −1.5 dex are shown with a density plot. There are two colour codes for the interlopers: dwarf stars are displayed in deep purple and isochrone-filtered sources in light purple. Left panel: before filtering, 9 432 286 stars. Right panel: after, 2 880 338 stars (Sect. 3.3). The blue curve corresponds to the dwarf and giant separation of Fernández-Alvar et al. (2021).</figcaption></figure>

<figure><img src="image_2.png" /><figcaption>**Fig. 2.** Spatial density distribution of the isochrone-filtered PGS in logarithmic scale (only giant stars). Left panel: Galactocentric Y versus Galactocentric X. Middle panel: Galactocentric Z versus Galactocentric cylindrical radius R. Right panel: $V_\phi$ versus [Fe/H]; the different ellipses qualitatively indicate the location of several identified Galactic components, following Fig. 8 of Fernández-Alvar et al. (2021).</figcaption></figure>

# 4. Validations using spectroscopic surveys

## 4.1. Validation of the isochrone filtering procedure

To confirm the reliability of PGS, and to evaluate the amount of metal-rich contaminants, we used the publicly available spectroscopic metallicity catalogues of APOGEE DR17 (Majewski et al. 2017; Abdurro'uf et al. 2022), GALAH DR3 (De Silva et al. 2015; Buder et al. 2021), LAMOST DR8 (Cui et al. 2012; Chen et al. 2022), and Gaia RVS GSP-spec (Recio-Blanco et al. 2016, 2023). It is important to note that in this section, we kept both dwarfs and giants, although white dwarfs were discarded for all surveys (stars with $M_{G_0} \ge 7.5$). The orbital parameters of the spectroscopic catalogues were computed similarly to the method described in Sect. 2. For APOGEE and GALAH, we used their spectroscopic RVs as input, whereas we used Gaia RVs for LAMOST and GSP-spec. Table B1 gives the sample size of each cross-match before and after performing the isochrone filtering method. The following

quality cuts were applied to each cross-match between PGS and the spectroscopic catalogues:

- APOGEE: following Abdurro'uf et al. (2022) and Hegedűs et al. (2023), we removed stars with any of its parameter flagged as _BAD<sup>11</sup>;

- GALAH: we followed Buder et al. (2021) and Hegedűs et al. (2023) and the best practices for the use of GALAH DR3<sup>12</sup>, that is, columns and flags recommendations;

- LAMOST: following Luo et al. (2015) and Tsantaki et al. (2022), we selected stars with $S/N > 15$ in the g- and i-bands at all [Fe/H] ranges (see Li et al. 2018 for limitations), and excluded stars with negative RV errors;

- GSP-spec: following Table 2 of Recio-Blanco et al. (2023), we selected stars with null values for flags 1–13 (related to

<footnote n="11"><sup>11</sup> https://www.sdss4.org/dr17/irspec/apogee-bitmasks/#APOGEEBitmasks</footnote>

<footnote n="12"><sup>12</sup> https://www.galah-survey.org/dr3/using_the_data/</footnote>
atmospheric parameters except abundances)<fnref n="13" />, except flags 7 (fluxNoise) and 8 (extrapol) for which values ≤2 were allowed. Table C1 summarises the different flags applied to this sample. Additionally, some hot metal-poor stars were discarded as they are likely badly parameterised (i.e., with $T_{\text{eff}} > 6000$ K and $[\text{Fe/H}]_{\text{GSP-spec}} < -1$, see Sect. 10.5 of Recio-Blanco et al. 2023). We note that for GSP-spec, we use $[\text{M/H}]_{\text{calibrated}}$ as our metallicity estimator.

In the following, we focus on comparisons with APOGEE since it is of high resolution and highly reliable for giant stars. Besides, APOGEE contains a large number of higher-metallicity disc stars, our primary source of contaminants. Validation using all spectroscopic catalogues can be found in Appendix B.

The left panel of Fig. 3 shows the cross-match between PGS and APOGEE in metallicity, where the broad majority of the sample follows a one-to-one relation. Between -2.5 < [Fe/H]_{APOGEE} < -1, we see a relatively good agreement. APOGEE metallicities do not go below -2.5 because of the lower boundary of their adopted training grid (García Pérez et al. 2016), whereas some PGS stars can be as metal-poor as -4.0. Thus, for these stars, APOGEE values are most likely biased, while PGS values can still be reliable. However, a fraction of PGS metal-poor candidates ([Fe/H]_{PGS} < -1.5) have nearly solar APOGEE metallicities. This is manifesting as a horizontal sequence in the left panel of Fig. 3 between -4.0 < [Fe/H]_{PGS} < -1.5, and makes up for 8.04% of the total number of PGS metal-poor candidates within the PGS-APOGEE cross-match. In what follows, we dub these contaminants ‘interlopers’. Although in small numbers, the metal-rich interlopers evidenced most likely belong to the disc and might contaminate highly prograde VMP stars. For this paper that aims to characterise VMP disc stars, the number of interlopers therefore needs to be minimised, which is why they are removed from our final sample.

To infer whether a certain type of star causes the bias observed between both metallicity estimates, we looked at the CMD of the PGS-APOGEE cross-match (left panel of Fig. 4). The majority of the metal-poor subsample is located on the blue side of the RGB, while interlopers can be found either on the main sequence or on the red side of the RGB. We note that interlopers are unlikely to be APOGEE stars with a wrong colour estimation (therefore, a wrong $T_{\text{eff}}$), as no clear dependency on the extinction coefficient $E(B-V)$ was identified. The isochrones with different combinations of age and metallicity visually translate the arguments used to construct the isochrone filtering method (Sect. 3.3). After applying the isochrone filtering method, 94% of APOGEE interlopers are removed (see middle panel of Fig. 3 and right panel of Fig. 4). For LAMOST, GALAH and GSP-spec, the percentages are detailed in Appendix B. The strength of this method can be seen in the middle panel of Fig. 3, where the majority of the horizontal sequence caused by interlopers was successfully removed.

We note that this method does prevent some interlopers from being discarded because they occupy the region where metal-poor stars physically exist, and conversely, it removes badly parameterised stars that are not interlopers but that still occupy the region where metal-poor stars should not exist. The latter represents an excess of 2.4% of the PGS metal-poor sample that might have been removed. Finally, it is also possible that outliers from APOGEE make up for a fraction of the undiscarded percentage of interlopers. In the right panel of Fig. 3, the aftermath of removing dwarf stars from our sample using the separation criterion of Fernández-Alvar et al. (2021) is the disappearance

of the overdensity at $[\text{Fe/H}]_{\text{PGS}} \sim -0.3$, and an overall lower density along the [1:1] diagonal.

## 4.2. Estimating contamination in the isochrone filtered PGS

After applying the isochrone filtering method to our PGS sample, we statistically assessed the likely remaining contamination using scaled comparisons with the spectroscopic catalogues. Here, we investigated the PGS-APOGEE cross-correlation. Other surveys are discussed in Appendix D. We adopted the following procedure:

- for a given $[\text{Fe/H}]_{\text{PGS}}$ range, we counted and flagged as ‘spectroscopic contaminant’ every star (already isochrone-filtered) located beyond the $2\sigma$ threshold of the $\Delta_{[\text{Fe/H}]} = [\text{Fe/H}]_{\text{PGS}} - [\text{Fe/H}]_{\text{APOGEE}}$ distribution (see Fig. 5);

- we then rescaled the contaminants: for a given $[\text{Fe/H}]_{\text{PGS}}$ range, we multiplied the number of contaminant counts by a factor, corresponding to the total number of PGS stars in the given $[\text{Fe/H}]_{\text{PGS}}$ range, divided by the number of PGS stars with spectroscopic estimates (including the contaminants), in the same $[\text{Fe/H}]_{\text{PGS}}$ range.

This definition of contaminants is conservative because it includes a significant number of true metal-poor stars. Specifically, it may include stars with underestimated spectroscopic metallicities as well as stars with only slightly overestimated photometric metallicities, which still belong to the metal-poor regime and cannot be mistaken for metal-rich disc contamination.

## 4.3. Impact of the spectroscopic contamination

In this section, we analyse the $V_\phi$ distribution and the action space of the isochrone-filtered PGS after estimating the spectroscopic contamination using APOGEE.

### 4.3.1. In the $V_\phi$ distribution

Figure 6 shows the $V_\phi$ distribution of PGS with associated spectroscopic contaminants in three $[\text{Fe/H}]$ intervals. In the top panel, for the most metal-rich interval ($-1.0 < [\text{Fe/H}] < -0.3$), the bulk of the PGS sample (in cyan) is compatible with the circular velocity measurement in the solar vicinity (mostly thin disc stars) of ~230 km.s$^{-1}$ (Eilers et al. 2019; Vieira et al. 2022; Pöder et al. 2023). The population is mildly affected by contaminants (in purple). The middle panel, with $-1.7 < [\text{Fe/H}] < -1.0$, is interesting because the PGS distribution seems to be compatible with three components: one centred at ~0–10 km.s$^{-1}$, in line with a stationary and, or slightly prograde halo (Tian et al. 2019), a second one at ~180 km.s$^{-1}$, compatible with the thick disc value of Fernández-Alvar et al. (2021), and a third one centred at ~240–250 km.s$^{-1}$, similar to the one identified in the top panel and dominated by metal-rich contamination, as expected at high $V_\phi$ (Hayden et al. 2015; Kordopatis et al. 2015). Finally, in the bottom panel, where the most metal-poor stars reside ($-4.0 < [\text{Fe/H}] < -1.7$), the PGS sample is centred at 0 km.s$^{-1}$, although a prograde-retrograde asymmetry is definitely visible, with a skewed distribution towards prograde values. Contaminants are spread homogeneously along the entire $V_\phi$ range, and are not discarding PGS VMP stars with high $V_\phi$. Out of all the metallicity intervals, the one between -1.0 and -1.7 is the most contaminated. The reason for the presence of higher contamination rates at high $V_\phi$ in this particular interval can be seen in Fig. 5. The middle panel displays a strong asymmetry in the

<footnote n="13">13 flags_gspspec LIKE '00000000000000%'.</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 3.** APOGEE versus PGS metallicities, colour-coded in number density. Left panel: before applying the isochrone filtering method (dwarfs + giants). Middle panel: after applying the isochrone filtering method (giants + dwarfs). Right panel: after applying the dwarf and giant separation criterion of Fernández-Alvar et al. (2021) (only giants, see Sect. 4.1). The black dashed line represents the [1:1] relation.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 4.** CMD of the PGS-APOGEE cross-match (grey). The metal-poor subsample ([Fe/H]$_{PGS}$ < −1.5 dex), is colour-coded in cyan. Interlopers are colour-coded in purple (Sect. 4.1). PARSEC isochrones of different age-metallicity combinations are also displayed on the right panel.</figcaption></figure>

distribution of contaminants with respect to the two side panels. In fact, this contamination is mostly made of metal-rich (solar) stars shifted in [Fe/H] down to almost −1.5; as explained in Sect. 1, this is because the metal-rich population is dominant in the thin disc and is more likely to pollute VMP prograde samples. Therefore, it is not surprising to find contamination at MW thin disc $V_\phi$.

### 4.3.2. In the action space

Figure 7 shows the action space of our PGS sample in the same [Fe/H] intervals chosen in Fig. 6. It is colour-coded in density, and contour lines of the 33, 66 and 90% of the sample are also displayed.

The action space is defined as $(J_z - J_r)/J_{tot}$ versus $J_\phi/J_{tot}$, where $J_r$ is the radial action, $J_\phi$ is the azimuthal action ($J_\phi = L_z$ since we assume an axisymmetric potential), $J_z$ is the vertical action, and $J_{tot}$ is $J_r + |L_z| + J_z$. The action space allows us to identify different families of orbits (Binney 2012a,b):

- radial: $L_z/J_{tot} = 0$ and $(J_z - J_r)/J_{tot} = -1$; the radial component dominates

- polar: $L_z/J_{tot} = 0$ and $(J_z - J_r)/J_{tot} = +1$; the vertical component dominates

- prograde: $L_z/J_{tot} \to 1$

- retrograde: $L_z/J_{tot} \to -1$

- circular: $J_z/J_{tot} > J_r/J_{tot}$

- planar: $J_z/J_{tot} < J_r/J_{tot}$

In the left panel of Fig. 7, 66% of the sample is labelled as circular-prograde-planar, which undoubtedly corresponds to a thin disc population with high $V_\phi$ (hence high angular momenta $L_z$), low eccentricities (circular orbits) and low $Z_{max}$ (close to the plane). More specifically, the mean eccentricity of the 66% of the sample is 0.20 and the mean $Z_{max}$ is 1.3 kpc. The rest of the sample extends along the planar down to radial regions of the action space with null $L_z$; this structure could come from different contributions, mainly the *Splash* (Belokurov et al. 2020) or the thick disc. In the middle panel, 33% of the PGS sample is prograde-planar between −1.7 < [Fe/H] < −1.0, and can be associated with the thin disc. However, it is difficult to completely trust this distribution in the action space, due to the high contamination rates at high $V_\phi$ (middle panel of Fig. 6). In the right panel, the bulk of the sample extends towards the planar and radial regions. The 90% contour line expands more towards retrograde values than for the two previous metallicity ranges; this is in line with the second peak in the $V_\phi$ distribution that we attribute to a halo contribution. The prograde-planar population is found again in the right panel, although the retrograde and prograde regions are more evenly distributed. Altogether, the error bars in the action space are small with respect to the signatures that we seek, especially since they are even smaller for prograde stars. Errors are driven mostly by distance errors, which stem from the maximum errors on parallaxes allowed in our selection. We note that the depletion of data when applying a stricter cut on the parallax (e.g. $\delta\varpi/\varpi < 0.10$) is localised in the polar area of the action space; however, we recovered the same trends as the study conducted with the original cut of 0.20.

To get a better understanding of what to expect in terms of dynamics, we compared our sample with the Gaia Universe Model Snapshot (GUMS, Robin et al. 2012), a model based on the Besançon Galaxy Model (BGM, Robin et al. 2003). This model generates a synthetic catalogue at a given static time (snapshot) to simulate the environment in which Gaia objects are or will be observed. It is semi-empirical as it is partly based on stellar evolution theory and constrained by observations. GUMS contains the main Galactic components (the thick and thin discs,
<figure><img src="image_5.png" /><figcaption>**Fig. 5.** $\Delta_{[Fe/H]} = [Fe/H]_{PGS} - [Fe/H]_{APOGEE}$ distribution (black). The Gaussian fit to the distribution is overplotted in red dashed lines. $\sigma$ and $2\sigma$ limits are respectively indicated in blue and teal vertical lines. The PGS-APOGEE contaminants, i.e. sources with $\Delta_{[Fe/H]} = [Fe/H]_{PGS} - [Fe/H]_{APOGEE} > 2\sigma$ are colour-coded in purple. Left panel: $-1.0 < [Fe/H]_{PGS} < -0.3$. Middle panel: $-1.7 < [Fe/H]_{PGS} < -1.0$. Right panel: $-4.0 < [Fe/H]_{PGS} < -1.7$ (Sect. 4.2).</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Fig. 6.** $V_{\phi}$ distribution of PGS. Counts with associated error bars (Poisson noise) are displayed in blue while PGS-APOGEE spectroscopic contaminants are displayed in purple. Contaminants are scaled by a factor equal to the size of PGS in each $[Fe/H]_{PGS}$ interval, divided by the size of APOGEE in the same intervals. Top panel: $-1.0 < [Fe/H]_{PGS} < -0.3$. Middle panel: $-1.7 < [Fe/H]_{PGS} < -1.0$. Bottom panel: $-4.0 < [Fe/H]_{PGS} < -1.7$.</figcaption></figure>

the halo, the bulge and the bar), but no stellar streams, nor any other known substructures. The parameters chosen for the Galactic components are constrained by the BGM<fnref n="14" />. For this mock model, we considered GUMS up to $G = 15.3$ mag, for $|b| > 7$ deg, and for giant stars (defined here as $T_{\text{eff}} < 6000$ K and $\log g < 3.9$ dex) to be able to compare with our catalogue of PGS giants. We also removed very young stars ($<0.15$ Gyr) since, as indicated in the Gaia documentation, such stars have had their kinematics badly implemented into the simulated catalogue. We computed stellar orbits using the McMillan (2017) potential, in the same fashion as in Sect. 2; we note that, to remain coherent with the GUMS dynamical model, we adopted the GUMS solar velocities (12.75, 0.93, 7.10) km.s⁻¹ (Robin et al. 2017), as well as $V_{\text{LSR}} = 230.6$ km.s⁻¹ at $R_{\odot} = 8$ kpc (Sofue 2015). In this

first selection function of GUMS stars, we verified that the main Galactic components were retrieved as expected in the BGM. The thin disc contains stars with a maximum scale length of 2011 pc and a mean rotation of 221.85 km.s⁻¹. The thick disc is modelled as a ‘young’ and ‘old’ thick disc of respectively 10 and 12 Gyr, in respective proportions of 68 and 32%; we found the expected scale heights of 400 and 795 pc, and scale lengths of 2040 and 2919 pc, and a mean rotation for the entire population 186.52 km.s⁻¹. Finally, halo stars have an eccentricity of 0.774, and a mean rotation of −0.34 km.s⁻¹. We stress that although GUMS does not capture the fine details of Galactic structure that would be needed to properly model our metal-poor population, it is a very helpful guide to understand the behaviour of classical MW stellar populations in various kinematical spaces, since it has been shown to reproduce the overall star counts and kinematics of the Galaxy.

We display the GUMS action space in Fig. 8. Out of the sample of GUMS thin disc stars, 90% are located on highly prograde and circular orbits; the GUMS thick disc is essentially prograde with equal amounts of radial and polar orbits. The GUMS halo is very symmetrical with respect to $L_z/J_{\text{tot}} = 0$. As expected, the left panel of Fig. 8 resembles the distribution of our data in the left panel of Fig. 7. According to the distribution of the GUMS thin and thick disc, we suspect that our sample at intermediate $[Fe/H]$ is a mix of thin and thick disc stars. As opposed to the symmetry observed for the GUMS halo, our sample at the lowest $[Fe/H]$ interval (right panel of Fig. 7) is clearly asymmetric in favour of the prograde component. Given the high contamination of PGS at intermediate $[Fe/H]$ (between −1.7 and −1.0) by disc stars, it would be risky to interpret the high $V_{\phi}$ PGS VMP stars ($[Fe/H] < -1.7$) without considering contamination. Although it is low overall in that range (see right panel of Fig. 5), it may contain a higher proportion of stars with disc kinematics than the truly metal-poor parent sample. In the following section, we statistically decontaminate our PGS sample using spectroscopic estimates.

## 5. Statistical decontamination of PGS

In this section, we detail how we decontaminated PGS for $-4.0 < [Fe/H] < -1.7$, using spectroscopic estimates of contamination. We worked with APOGEE and LAMOST in the $V_{\phi} - [Fe/H]$ space and in the action space. As discussed in Appendix D, due to completeness (GSP-spec) and statistical (GALAH) issues, we investigated the decontamination with APOGEE because of its high resolution, and with LAMOST because of its statistics, despite its low resolution. To remain within the range of 0.3 dex in $[Fe/H]$ uncertainties, we computed the mean of contamination

<footnote n="14">¹⁴ https://gea.esac.esa.int/archive/documentation/GEDR3/Data_processing/chap_simulated/sec CU2UM/ssec CU2starsgal.html</footnote>
<figure><img src="image_7.png" /><figcaption>**Fig. 7.** Action space: ($J_z - J_r$) as a function of $L_z$. Axes are normalised by $J_{tot} = |L_z| + J_r + J_z$. Left panel: $-1.0 < [\text{Fe/H}]_{\text{PGS}} < -0.3$. Middle panel: $-1.7 < [\text{Fe/H}]_{\text{PGS}} < -1.0$. Right panel: $-4.0 < [\text{Fe/H}]_{\text{PGS}} < -1.7$. We note that the colour bar range changes for each interval. We add mean error bars for each quadrant of the action space in light pink.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Fig. 8.** Gaia Universe Model Snapshot (GUMS) action space. Left panel: thin disc stars. Middle panel: thick disc stars. Right panel: halo stars. The contour lines at 33, 66, 90 and 99% are plotted in black dashed lines. We note that the colour bar range changes for each panel.</figcaption></figure>

in different areas of the $V_\phi$ - [Fe/H] space and of the action space. In each area, the scaling factor for contamination changes with the density of stars in PGS, and the density of stars in the survey, in designated areas of the $V_\phi$ - [Fe/H] space and of the action space. For both surveys, we checked :

- the density of stars in the contaminants samples;

- PGS minus individual contamination (pixel-by-pixel);

- PGS minus the mean of contamination.

## 5.1. $V_\phi$ - [Fe/H] space

Fig. 9 summarises the main steps of decontamination in the $V_\phi$ - [Fe/H] space. Here, we plotted 30 individual areas, ranging from -250 to 250 km.s⁻¹ vertically (spaced by 50 km.s⁻¹), and from -1.7 to -4.0 horizontally (spaced at -2.0 and -2.5). We note that the vertical bin size is 10 km.s⁻¹ and the horizontal bin size is 0.1 dex. The first panel (top) shows the density of the original PGS data. We notice an asymmetry in $V_\phi$ which was already observed in the bottom panel of Fig. 6. For [Fe/H] < -2.5, there is a sharp decrease in density, reflecting the shape of the MDF of the MW; even with these low counts, an asymmetry in the $V_\phi$ distribution is visible. The second row of panels shows the density of the spectroscopic contamination. Starting with the left panel (APOGEE), we notice blank areas at [200, 250] and [-150, -250] km.s⁻¹; we have no information on the con-

tamination there since APOGEE does not provide spectroscopic metallicity measurements below -2.5. Hence, only contaminants with [Fe/H]$_{PGS}$ < -2.5 and [Fe/H]$_{APOGEE}$ > -2.5 can appear in that range. When there are none, no contamination can be measured. Nonetheless, for the stars that verify this condition, we have evidence of a low contamination rate (at [150, 200] and [0, -50] km.s⁻¹). Overall, the contamination is lower with decreasing [Fe/H], but it is also very high in the three upper right areas, corresponding to the population we are aiming to characterise. In the right panel, as expected, we are able to follow the contamination down to -4.0 with LAMOST because of the more complete sample available. However, this time contamination does not decrease with decreasing [Fe/H], but is rather concentrated between -1.7 and -2.5 dex. It is not restricted to high $V_\phi$ values but also extends to the [0, -50] km.s⁻¹ area. The third row of panels shows the pixel-by-pixel decontamination of PGS. We note that the blank APOGEE contamination areas were replaced by PGS densities, since there is in practice no contamination to account for. In both cases, the decontamination is causing pixels with negative values, due to statistical fluctuations of the low contamination rates (the contaminant pixel removed is too dense in comparison with the PGS pixel). This is highly dependent on the sample sizes of PGS, APOGEE and LAMOST since they rule the weight of the scaling factor. However, this view is useful to identify the areas most affected by the decontamination. In particular, we note that regardless of the negative structure, the
<figure><img src="image_9.png" /><figcaption>**Fig. 9.** Densities of stars in the $V_\phi$ – [Fe/H] space. The grid corresponds to areas delimited in $V_\phi$ and [Fe/H]. The scaling of contaminants varies from one area to another with the density of PGS in each area. The dimension of a pixel is $0.1 \text{ dex} \times 10 \text{ km.s}^{-1}$. First panel (top): PGS only. In the next rows of panels: left = APOGEE, right = LAMOST. Second row: density of spectroscopic contaminants, scaled to the sample size of PGS in each pixel. Third row: PGS decontaminated pixel-by-pixel, i.e. subtraction between a pixel of PGS density and a pixel of spectroscopic contamination density. Fourth row: residuals percentages of the pixel-by-pixel decontamination of PGS. Fifth row: PGS decontaminated by the mean of the spectroscopic contamination, i.e. by the mean of the spectroscopic contamination density in a designated $V_\phi$ – [Fe/H] area.</figcaption></figure>

density of PGS in the upper right areas of the $V_\phi$ – [Fe/H] space ($V_\phi > 150 \text{ km.s}^{-1}$, $-2.5 < [\text{Fe/H}] < -1.7$) is still remarkably high when decontaminating with APOGEE, and non-negligible when decontaminating with LAMOST.

Another way to visualise the effect of the pixel-by-pixel decontamination is to show the residuals, after removing spectroscopic contaminants, as seen in the fourth row of panels. After

both decontaminations, we notice that between 50 and 100% of the sample is conserved, especially at $[-100; 100] \text{ km.s}^{-1}$ for all metallicities. We do see more clearly in the upper right areas of the $V_\phi$ – [Fe/H] space mentioned before that a large percentage of data is preserved.

To overcome the statistical fluctuations, we also display the result in the fifth row of panels, when the decontamination is
made by subtracting the mean value of contamination in each area, delimited the dashed lines corresponding to different combinations of $V_φ$ and [Fe/H]. Impressively, there is not much difference between both decontaminations, and also between the decontaminations and the original PGS. Again, we do see a change in the densities of the upper right areas, but the population of VMP planar-prograde stars we suspected could be erased after decontaminating is in fact still present.

## 5.2. Action space

The principle of decontamination in the action space is the same as in the $V_φ$–[Fe/H] space. The main difference is that we manipulated the action space to resemble the layout of Fig. 9, so that comparisons between both spaces are more straightforward. To do so, we rotated the (central) coordinates of each pixel composing the original action space (top panel of Fig. 10) by +45 degrees. Then, we assigned a weight to each pixel, equal to the scaling factor defined in Sect. 4.2. The limits of this newly defined action space correspond to the rotated coordinates of each corner of the original action space; the resulting space is displayed in the bottom panel of Fig. 10. The rotated action space is subsequently divided into nine even areas of length ~0.4 in both directions.

Figure 11 summarises the decontamination in the action space with APOGEE. The results of the decontamination with LAMOST, which are very similar to APOGEE, are displayed in Fig. E1. Comparably to Fig. 9, regardless of the method employed and the survey used to decontaminate PGS, we notice a high fraction of stars with circular-prograde-planar orbits, and also prograde-planar orbits in the decontaminated PGS.

# 6. VMP stars on disc-like orbits

After decontaminating PGS in the $V_φ$ – [Fe/H] space and in the action space, regardless of the survey used there are strong hints indicating the presence of a population with the following characteristics:

- high $V_φ$ ($>180$ km.s$^{-1}$);

- prograde-planar orbits (10% of the PGS VMP sample in the rotated action space, between 0.2 and 0.7, regardless of the survey used for decontamination);

- metallicities below −1.7 dex.

In this section, we aim to provide a kinematical and dynamical characterisation of the highly prograde<fnref n="15" /> and planar population ($Z_{max} < 1.5$ kpc), which makes up for 2% of the PGS VMP sample, and to constrain its lower limit in [Fe/H], with special attention to test whether this is a tail of the MDF of known stellar populations (thin and thick disc, prograde halo).

## 6.1. Trends in the spatial distribution

### 6.1.1. X-Y and R-Z

First, we investigated the spatial distribution of the isochrone-filtered PGS<fnref n="16" /> in Fig. 12. As previously noted in Sec. 3.3 with Fig. 2, there is a visual asymmetry in the X-Y and R-Z distribution of PGS, where stars are mostly located in the inner parts

<figure><img src="image_10.png" /><figcaption>**Fig. 10.** Densities in the action space. Top: original PGS after applying the isochrone filtering method, in the most metal-poor interval (−4.0 ≤ [Fe/H] < −1.7). Same plot as the right panel of Fig. 7. Bottom: same as the top panel, after rotating the (central) coordinates of each pixel by +45 degrees. The limits of the plot in that space correspond to the rotated coordinates of the four corners of the original action space. The grey dashed lines of the bottom panel correspond to different values of $L_z/J_{tot}$: [−0.95, −0.8, −0.5, 0, 0.5, 0.8, 0.95].</figcaption></figure>

of the Galaxy. In Fig. 12, we show the full PGS with a selection of disc-like stars, that is, stars with $V_φ$ > 180 km.s$^{-1}$ and $Z_{max} < 1.5$ kpc (density plot). We also plot a subsample of disc-like stars with [Fe/H] < −1.7 dex (teal scatter). We note that, as expected from our cut in $Z_{max}$, the R-Z distribution of the disc-like population is concentrated very close to the plane.

We can divide the R(X)-Z plane (right panel) into four quadrants to evidence any asymmetry in the full PGS sample, the disc-like subset and the VMP disc-like subset. First, we looked for asymmetries in Y (left panel); overall, the fraction of stars with $Y > 0$ kpc is relatively constant for the full PGS sample, its disc-like subset (both 43%), and the disc-like subset (38%) for which we evidenced a non-negligible asymmetry with negative Y stars.

In the R-Z plane, we did not evidence strong asymmetries in R for the full PGS sample (52% of the sample lies below $R = 8$ kpc) nor Z (45% above 0 kpc). The disc-like subset also does not show strong asymmetries in R (47% of the subset below 8 kpc) nor Z (also 45% above 0 kpc). We conclude that the asymmetries identified in Fig. 2 are simply probable visual effects from the Gaia scanning law. We notice a pronounced asymmetry in R and Z for the VMP disc-like subset, with 65% located below $R = 8$ kpc, and 42% with positive Z. When checking R and Z simultaneously, 29% of disc-like stars have $R < 8$ kpc and $Z < 0$ kpc, while for the VMP disc-like subsample it represents

<footnote n="15">¹⁵ In the following, we refer to highly prograde stars as stars with $V_φ$ > 180 km.s$^{-1}$.</footnote>

<footnote n="16">¹⁶ In Sects. 6.1.1 and 6.1.2, we neglected the contamination of the PGS metal-poor sample, as we have shown that for metallicities below −1.7, it is not a dominant effect on our sample, and because we do not expect contamination to impact directly the spatial and $Z_{max} - R_{max}$ distributions examined in these sections.</footnote>
<figure><img src="image_11.png" /><figcaption>**Fig. 11.** Densities in the action space, in the most metal-poor interval (−4.0 ≤ [Fe/H] < −1.7). The grid corresponds to areas in the action space delimited in ($J_z - J_r$) and $L_z$. From top to bottom: first panel: density of APOGEE contaminants. Second panel: PGS decontaminated pixel-by-pixel, using APOGEE. Third panel: residuals percentages of the pixel-by-pixel PGS decontamination using APOGEE. Fourth panel: PGS decontaminated by the mean of the APOGEE contamination. The grey dashed lines correspond to different values of $L_z/J_{tot}$: [−0.95, −0.8, −0.5, 0, 0.5, 0.8, 0.95].</figcaption></figure>

41%, a value consistent with Thomas et al. (2019). These observations could be hints of a disc population with a shorter scale length than the MW thin disc, thus reminiscent of the thick disc (Bensby et al. 2011; Bovy et al. 2012; Hayden et al. 2015), or the flaring of the disc (Thomas et al. 2019). It could also correspond to halo stars, for which the density is higher in the inner Galaxy (Rix et al. 2022).

### 6.1.2. $Z_{max} - R_{max}$

Following Haywood et al. (2018), Di Matteo et al. (2020), and Koppelman et al. (2021), we then probed the $Z_{max} - R_{max}$ plane. Here, $R_{max}$ corresponds to the square root of the difference between the apocenter cylindrical radius, $R_{apo}$, and $Z_{max}$. $R_{max} = \sqrt{R_{apo}^2 - Z_{max}^2}$. We used this space to better constrain the bulk of the distribution of each GUMS Galactic component (thin disc, thick disc, halo), in order to compare it with PGS, below [Fe/H] = −1.7. To achieve this, we selected for each PGS star, the nearest neighbour from GUMS in longitude $\ell$, latitude $b$, $\varpi$ and apparent G magnitude simultaneously. We allowed the procedure to select several times the same GUMS star if it is the closest match. We obtained:

- 13% of nearest neighbour repeats for the GUMS halo sample;

- 12% of repeats for the GUMS thick disc sample;

- 53% of repeats for the GUMS thin disc sample.

Interestingly, the amount of repeats in each component reflects the difference between the spatial distribution of our PGS metal-poor sample, and that of a pure halo, a pure thin or thick disc. In fact, the observed spatial distribution of our PGS metal-poor sample does not resemble a classical thin disc, as seen in Fig. 13; thus, forcing a matching selection in $\ell$, $b$, $\varpi$ and G magnitude results in a larger fraction of duplicates. In the following, we only worked with unique stars from our GUMS samples, therefore the GUMS halo, thick and thin disc samples are respectively 13, 12 and 53% smaller than our PGS VMP sample.

Figure 13 compares the distribution of PGS stars below −1.7 dex with the three GUMS Galactic components convolved by the selection function as described above. For legibility, the 66% contour lines of GUMS are plotted in red solid lines, and the 66% contour line of the PGS sample is plotted in yellow solid lines. We also added the respective contour lines of each highly prograde subsample ($V_\phi > 180 \text{ km.s}^{-1}$) in red and yellow dashed lines. The orange solid horizontal line corresponds to a constant $Z_{max} = 2 \text{ kpc}$. To add more meaning to this representation, Table 2 shows the fraction of stars located below $Z_{max} = 2 \text{ kpc}$, for each GUMS component as well as PGS (in different [Fe/H] intervals) for the full samples and their highly prograde subsets.

The left panel displays the density of the GUMS thin disc. As indicated by the orange solid line defined as $Z_{max} = 2 \text{ kpc}$, the bulk of the distribution is located below $2 \text{ kpc}^{17}$, whether we check the full sample (84%) or solely highly prograde stars (89%). The overall $Z_{max} - R_{max}$ distribution of PGS shows that the main 66% of the sample is constrained on scale heights $Z_{max} < 7 \text{ kpc}$. Similarly, $R_{max}$ is contained within 2 and 13 kpc. If we consider the full PGS sample, 22% fall below $Z_{max} = 2 \text{ kpc}$. As for the highly prograde subsample, this fraction increases to

<footnote n="17">17 $Z_{max} = 2 \text{ kpc}$ may seem high for a classical thin disc component such as the one assumed in GUMS; we verified that this is largely due to applying the PGS VMP sample selection function, preferentially selecting stars at low R and high Z, which in turn enter together with the orbits of the stars, to produce $Z_{max}$ and $R_{max}$, eventually affecting the distribution in that space.</footnote>
<figure><img src="image_12.png" /><figcaption>**Fig. 12.** Spatial distributions. Grey scatter plot: same as in Fig. 2. Density plot: disc-like stars ($V_{\phi} > 180 \text{ km.s}^{-1}$ and $Z_{\text{max}} < 1.5 \text{ kpc}$). Teal scatter plot: disc-like VMP stars. Left panel: X versus Y. Right panel: R versus Z.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Fig. 13.** $Z_{\text{max}}$ as a function of $R_{\text{max}}$. The density plot corresponds to GUMS stars that follow the selection function of the PGS VMP sample (below $[\text{Fe/H}] = -1.7$) in $\ell, b, \varpi$ and apparent G magnitude. The red solid line corresponds to the 66% contour line of each GUMS sample, while the yellow solid line is the 66% contour line of PGS VMP stars ($[\text{Fe/H}] < -1.7$). The red dashed line corresponds to the 66% contour line of each highly prograde ($V_{\phi} > 180 \text{ km.s}^{-1}$) GUMS subsample; the yellow dashed line is the 66% contour line of PGS VMP highly prograde stars. The orange solid horizontal line corresponds to $Z_{\text{max}} = 2 \text{ kpc}$. Left panel: GUMS thin disc. Middle panel: GUMS thick disc. Right panel: GUMS halo.</figcaption></figure>

33% while the 66% contour line still rises to $Z_{\text{max}} \sim 6 \text{ kpc}$, showing that the PGS distribution is distinct (more extended in $Z_{\text{max}}$) from that of the GUMS thin disc.

The middle panel shows the density plot of the GUMS thick disc. The 66% GUMS contour lines are higher in $Z_{\text{max}}$ than in the previous panel. Interestingly, the 66% PGS contour line for all PGS VMP stars seems to match the outer contours (~90%) of the GUMS thick disc. We also notice that the contour line of the highly prograde PGS VMP subsample coincides with the contour line of the GUMS total thick disc. In fact, we found that 38% of GUMS thick disc stars lie below $Z_{\text{max}} = 2 \text{ kpc}$, which is only 5% more than the PGS VMP prograde subsample, a hint that the thick disc might significantly contribute to this population.

Finally, the right panel shows the density plot of the GUMS halo. By contrast to the previous panels, both 66% contour lines match with each other when considering the full samples. However, we note that the contour lines of the GUMS halo extend to higher values of $Z_{\text{max}}$, in particular for $R_{\text{max}} \le 8 \text{ kpc}$, where PGS VMP stars have significantly lower $Z_{\text{max}}$ than GUMS halo stars; furthermore, the highly prograde GUMS halo extends to $R_{\text{max}} > 15 \text{ kpc}$ (presumably including stars with very large $R_{\text{apo}}$), which is not observed in the PGS sample.

In summary, within the limitation of GUMS (i.e. a mock model with no additional substructures other than a prograde or dual halo, a thick disc and a thin disc), when no selection on $V_{\phi}$ is made, the comparison between GUMS components and

PGS leads to the conclusion that the best match in fraction is the halo. Simulating self-consistently a prograde halo (or a dual halo with one prograde component) is beyond the scope of this paper, although we have verified that a simple shift of $V_{\phi}$ in GUMS does not modify the $Z_{\text{max}} - R_{\text{max}}$ distribution. We therefore confirm that the majority of the PGS sample would be composed of stars with a halo $Z_{\text{max}} - R_{\text{max}}$ distribution, and halo kinematics, as we show in the following sections. However, when focusing only on the highly prograde subsample, the fraction of PGS stars below $Z_{\text{max}} = 2 \text{ kpc}$ regardless of the $[\text{Fe/H}]$ interval, is twice higher than the corresponding GUMS halo fraction, and resembles more the fraction seen in the GUMS thick disc. This is a hint that the spatial distribution and orbital configuration of the highly prograde PGS metal-poor population is consistent with that of a kinematically warm Galactic disc (possibly even slightly warmer than the canonical GUMS thick disc). In the next section, we investigate the $V_{\phi}$ distribution and the action space of the most metal-poor PGS stars, after decontaminating with spectroscopic surveys, to further characterise the kinematics and orbital configuration of this population.

## 6.2. Asymmetries in the decontaminated PGS

### 6.2.1. $V_{\phi}$ distribution

To identify and quantify a potential contribution from the halo to our population of interest, an empirical test is to fold the $V_{\phi}$
distribution after decontaminating PGS. This enables to visualise any remaining prograde count after subtracting retrograde counts. But, folding along an axis of symmetry centred on 0 amounts to assuming that the halo has zero net rotation, which is still an open debate (Morrison et al. 1990; Deason et al. 2011, 2017; Tian et al. 2019; Fernández-Alvar et al. 2021).

Therefore, we investigated the folded $V_φ$ distribution below -1.7 dex with increasing halo velocity central value, with the aim of finding the limit for which the highly prograde population is entirely removed; this is summarised in Fig. 14. Because the results for LAMOST are very similar, we show them in Fig. E2. The first two (top) panels show the decontamination by APOGEE in the $V_φ$ distribution, similar to the pixel-by-pixel decontamination detailed in Sect. 5. The third panel shows the folding in the case where we assume a stationary halo, centred at $0 \text{ km.s}^{-1}$. The red dashed lines delimit the $V_φ$ at which we expect our population to be (180–300 km.s$^{-1}$). In both decontaminated and folded PGS $V_φ$ distributions, high $V_φ$ stars are present in non-negligible numbers (maximum counts ~250 for APOGEE and 200 for LAMOST). This means that removing the entire contribution from halo stars when the halo is stationary is not sufficient to suppress highly prograde VMP stars. We note that stars with thick disc-like $V_φ$ (~100 km.s$^{-1}$) are also present, in high numbers. In the following fourth to the eighth panels, we progressively shift the central value of the halo by $10 \text{ km.s}^{-1}$ until we reach $50 \text{ km.s}^{-1}$. The value of $30 \text{ km.s}^{-1}$ (panel 6) simultaenously minimises the global counts and does not over-correct from the prograde component (negative structure) at high $V_φ$ (180–300 km.s$^{-1}$), both for APOGEE and LAMOST. If we do not account for negative structure at high $V_φ$, the value of $40 \text{ km.s}^{-1}$ (panel 7) is the one that globally minimises the prograde population of our sample for both spectroscopic surveys<fnref n="18" />. Therefore, we conclude that a single-component halo would need to be prograde with $\bar{V}_φ$ of at least $30–40 \text{ km.s}^{-1}$ for the population of VMP highly prograde stars to be considered as the high-rotation tail of such a halo.

This value is in qualitative agreement with Zhang et al. (2023), who work with a sample of stars with Gaia RVs and Gaia DR3 XP spectra-based metallicities from Andrae et al. (2023). Their approach is to fit three-dimensional Gaussian Mixture Models (GMM) in a ($V_φ$, $V_R$, $V_Z$) space in bins of decreasing [Fe/H], from which they consistently find a dual halo; one that is stationary and one that rotates. If we combine the mixture and associated $V_φ$, the qualitative agreement becomes quantitative at $-2.0 < [\text{Fe/H}] < -1.6$ (24% of stationary halo rotating at $-8 \text{ km.s}^{-1}$ and 60% of prograde halo rotating at $\sim 72 \text{ km.s}^{-1}$ give an average rotation of $40 \text{ km.s}^{-1}$) and can also be verified below $-2.0$ dex, where only 28% of their population rotates. Taken individually, their prograde halo below $-1.6$ dex (two most metal-poor bins) rotates $\sim 30–40 \text{ km.s}^{-1}$ faster ($72–80 \text{ km.s}^{-1}$) than the one we infer here ($30–40 \text{ km.s}^{-1}$). This difference in mean rotation is most probably driven by the different approaches of our works. More specifically, our approach has been purely empirical; it is agnostic on the functional form (no assumption of gausssianity) and does not rely on any modelling of the velocity distributions. However, it is not free of hypothesis. First, we assumed that the highly retrograde and prograde tails of the velocity distribution are symmetrical; second, we assumed that the halo is a single component. This latter hypothesis is different from other works, including Belokurov & Kravtsov (2022) and Zhang et al. (2023) in a similar sample, or Ardern-Arentsen et al. (2024) in

<figure><img src="image_14.png" /><figcaption>**Fig. 14.** $V_φ$ distribution for stars with PGS metallicities between -1.7 and -4.0. From top to bottom: panels 1 and 2: PGS decontamination by APOGEE. Panels 3–8: folding of the retrograde counterparts over the prograde counterparts of the decontaminated PGS. The black dashed line indicates the location of our presumed halo. In panel 3, the folding is done assuming a theoretical halo. In the remaining panels, it is done assuming a prograde halo shifted by $10 \text{ km.s}^{-1}$ at every panel. The red dashed lines delimit the location of the population we aim to characterise.</figcaption></figure>

the inner Galaxy. The difference in the mean rotation of the halo needed to account for the population of highly prograde VMP stars may stem from this: a faster rotation of a less numerous sub-population (i.e. a second halo population) would be needed to erase that same signature.

### 6.2.2. Action space

To visualise how these results translate dynamically, we folded the rotated action space in Fig. 15. The trends obtained for LAMOST are displayed in Fig. E3; they are very close to the observations we make with APOGEE. The middle panel of Fig. 15

<footnote n="18">¹⁸ We verified that the results regarding the folding in the action space were not affected by a choice of either value.</footnote>
<table><caption>**Table 2.** Fraction of stars below $Z_{max} = 2$ kpc, for each GUMS component, and for PGS in different [Fe/H] intervals (all < −1.7 dex).</caption><thead><tr><th>Sample</th><th>Z<sub>max</sub> &lt; 2 kpc<br>out of all stars (%)</th><th>Z<sub>max</sub> &lt; 2 kpc out of the highly<br>prograde (V<sub>φ</sub> &gt; 180 km.s<sup>-1</sup>) sample (%)</th></tr></thead><tbody><tr><td>GUMS thin disc</td><td>84</td><td>89</td></tr><tr><td>GUMS thick disc</td><td>38</td><td>57</td></tr><tr><td>GUMS halo</td><td>18</td><td>14</td></tr><tr><td>PGS, [Fe/H] &lt; −1.7</td><td>22</td><td>33</td></tr><tr><td>PGS, −2.0 &lt; [Fe/H] &lt; −1.7</td><td>24</td><td>38</td></tr><tr><td>PGS, −2.5 &lt; [Fe/H] &lt; −2.0</td><td>21</td><td>30</td></tr><tr><td>PGS, −4.0 &lt; [Fe/H] &lt; −2.5</td><td>18</td><td>26</td></tr></tbody></table>

**Notes.** Each GUMS subset follows the selection function of the PGS VMP sample (below [Fe/H] = −1.7) in $ℓ, b, \sigma$ and apparent G magnitude. The fractions are given for full samples, and for highly prograde subsets ($V_φ > 180$ km.s$^{-1}$).

<figure><img src="image_15.png" /><figcaption>**Fig. 15.** Rotated action space between −4.0 ≲ [Fe/H] < −1.7. The grey dashed lines correspond to different values of $L_z/J_{tot}$: [−0.95, −0.8, −0.5, 0, 0.5, 0.8, 0.95]. Left column: decontamination process and folding using APOGEE contamination. Top row: density of spectroscopic contaminants. Middle row: pixel-by-pixel decontamination. Bottom row: folding of the retrograde components over the prograde components, for PGS decontaminated pixel-by-pixel. In this space, the folding axis of symmetry is the grey dashed line at $L_z/J_{tot} = 0$. Right column: same as the left panels, but the three actions were computed shifting the observed $V_φ$ by 40 km.s$^{-1}$.</figcaption></figure>

shows PGS decontaminated by APOGEE, pixel-by-pixel. The decontaminated PGS has a high density of stars with prograde-planar-circular orbits according to the grey dashed lines corresponding to $L_z/J_{tot} = 0.8$ and 0.95. But this structure is statistically significant along the prograde-planar region, down to the radial area. When folding the action space (prograde minus retrograde) with respect to the central diagonal (dashed line at $L_z/J_{tot} = 0$), we obtain the plot at the bottom left panel of Fig. 15, which is the equivalent of the third panel of Fig. 14. After folding, pixels with negative values appear below $L_z/J_{tot} = 0.5$. However, the prograde-planar region is still populated in relatively high counts. In particular, above $L_z/J_{tot} = 0.8$ the counts are the highest. This is in line with the observations of Fig. 14.

We now focus on the right column of Fig. 15. To model the effect of a shifted prograde halo in the action space, we

recomputed the three action components assuming a shifted $V_φ$, while $V_R$ and $V_Z$ remain unchanged. To be consistent with the limit value for the presence of VMP highly prograde stars found empirically, we shifted it by 40 km.s$^{-1}$. It is important to note that we folded the shifted retrograde counterparts onto the original prograde PGS counts. Since our main motivation is to get an idea of how the original action space evolves with a shifted halo component, we compared:

- the folding of the original action space, that is, the subtraction of the original retrograde component from the original prograde component;

- the folding of the original action space, but replacing the original retrograde component by the shifted (by 40 km.s$^{-1}$) retrograde component (the prograde component remains the same).
<figure><img src="image_16.png" /><figcaption>**Fig. 16.** Cumulative [Fe/H] distribution of the normalised counts of the decontaminated PGS, in different $V_\phi$ intervals. Left panel: $200 < |V_\phi| < 300$ km.s$^{-1}$. Middle panel: $100 < |V_\phi| < 200$ km.s$^{-1}$. Right panel: $0 < |V_\phi| < 100$ km.s$^{-1}$.</figcaption></figure>

That is, we subtracted the upper triangle matrix (with respect to the diagonal at $L_z/J_{tot} = 0$) representing the halo shifted by $40 \text{ km.s}^{-1}$ (middle right panel of Fig. 15) from the lower triangle matrix (with respect to the diagonal at $L_z/J_{tot} = 0$) representing the decontaminated PGS prograde counts (middle left panel of Fig. 15). Therefore, the result of the folding in the lower right panel should be treated very cautiously at $L_z/J_{tot} = 0$ because there is no physical meaning in matching the original PGS prograde counts with counts that are also prograde but were shifted to 0. As expected from the $V_\phi$ distribution of Fig. 14, the overall star counts in the residuals, in the bottom right panel of Fig. 15, are close to 0, which verifies our methodology. In spite of that, subtracting retrograde stars from our prograde sample has left a significant signature in the action space. While VMP stars with polar orbits are subtracted in excess (at all $L_z/J_{tot} < 0.8$ or even higher), stars with planar and circular orbits remain (where $L_z/J_{tot} > 0.95$ and stars with $L_z/J_{tot} > 0.8$ which are slightly polar). This indicates that although the rotational signature of these VMP stars with high $V_\phi$ can be explained by a prograde halo, their orbits are more radial and less polar than those of our empirical rotating halo. However, we might reach the limits of our method, since only the $V_\phi$ are arbitrarily shifted, which is not the case for the other velocities (E is not conserved anymore). From this dual view on the $V_\phi$ distribution and the action space, we conclude that if the halo is shifted by $40 \text{ km.s}^{-1}$, the remaining population is mostly planar-prograde, with $V_\phi$ of the order of $100–150 \text{ km.s}^{-1}$, which is compatible with a metal-poor thick disc (Kordopatis et al. 2013). Besides, we can highlight a clear asymmetry in the distribution of the prograde and retrograde stars at low [Fe/H], in favour of the prograde counterparts.

## 6.3. Insights on the [Fe/H] lower boundary

The lower limit in [Fe/H] of our sample can be constrained by comparing its MDF in different $V_\phi$ slices, which is simply another view of the bottom left panel of Fig. 9. We note that this particular panel takes into account PGS decontaminated by the mean of APOGEE contaminants. We plotted the cumulative distribution of metallicities below −1.7 dex in Fig. 16. The left panel shows high $|V_\phi|$ stars, the middle panel shows intermediate $|V_\phi|$ stars and the right panel shows stars with $|V_\phi|$ close to 0. Regardless of the panels, between −1.7 and −2.5 dex, there are globally as many highly prograde (>200 km.s$^{-1}$) stars

<figure><img src="image_17.png" /><figcaption>**Fig. 17.** [Fe/H] distribution of the decontaminated PGS, below −1.7 dex for different ranges of $V_\phi$.</figcaption></figure>

as highly retrograde (<−200 km.s$^{-1}$), if not more in the middle panel. Yet, below −2.5 dex the trend changes in favour of the prograde cumulative counts, which seem to dominate down to the lowest metallicities. This is especially visible for high $|V_\phi|$ stars (left panel). We draw similar conclusions when comparing the intermediate prograde and retrograde distributions.

To get a better sense of the prograde-retrograde asymmetry, we could check their ratios as a function of [Fe/H]. However, below −3 dex we over-decontaminate PGS (as discussed in Sect. 5.1); the resulting counts are negative, which is unsuited for determining these ratios. Therefore, we can only confidently check the prograde-retrograde ratio down to −3 at most. We compared three categories of prograde and retrograde populations: the high $|V_\phi|$ (>200 km.s$^{-1}$), the intermediate $|V_\phi|$ (~150 km.s$^{-1}$), and the small $|V_\phi|$ (~50 km.s$^{-1}$). We found:

- [0; +100] versus [0; −100] km.s$^{-1}$: the prograde-to-retrograde ratio is of the order of 1.5 in favour of the prograde;

- [+100; +200] versus [−100; −200] km.s$^{-1}$: the ratio is ~2.5 in favour of the prograde;

- [+200; +300] versus [−200; −300] km.s$^{-1}$: the ratio is ~2.6 in favour of the prograde.

If we solely check the original distribution (Fig. 17) below –2.5 dex, assuming a 2σ error bar based on Poisson noise, the
highly prograde population (>200 km.s⁻¹) is reaching 0 at –2.9 dex, that is, 1.3 dex lower than the lower limit of the disc population of Zhang et al. (2023). However, they do detect highly prograde stars below –2 dex, but they associate them to an underlying halo distribution rather than a separate population. Since their result is based on GMM, it is possible that our population of highly prograde and planar (Zmax < 1.5 kpc) stars is not detected below –1.6 dex, because it represents only 2% of our PGS VMP sample. Nevertheless, our result is in line with Fernández-Alvar et al. (2024), who also find this population below –1.5 dex.

## 7. Summary and conclusions

Our goal was to investigate and characterise the kinematics and orbital properties of VMP ([Fe/H] < –1.7) MW stars, with special emphasis on stars with disc-like orbits, namely prograde (in particular highly prograde stars, with V<sub>φ</sub> > 180 km.s<sup>-1</sup>) and planar (Z<sub>max</sub> < 1.5 kpc). To achieve this, we combined Pristine-Gaia synthetic photometry, Gaia astrometry and Gaia RVS radial velocities. We successfully filtered PGS from metal-rich outliers with a method based on isochrones in the CMD, and constructed a sample of ~3 M stars with –4 < [Fe/H] < 0. To further refine the quality of our data, we statistically decontaminated the PGS VMP subsample containing ~36 000 stars, using the comparison to the APOGEE and LAMOST spectroscopic surveys, through a criterion related to the difference between photometric and spectroscopic metallicity estimates. We found that contamination has minimal impact on the kinematical properties of the VMP sample. Globally, we recovered the classical halo trends in orbital Z<sub>max</sub> – R<sub>max</sub>, where the halo prediction of GUMS matches the bulk of our sample distribution, and in the action space, where a large fraction of our data (90% between –1.7 < [Fe/H] < –1.0 and 66% below –1.7) is distributed along the retrograde area.

Yet, similarly to Sestito et al. (2020), Bellazzini et al. (2024), and Viswanathan et al. (2024b) (see their Fig. 7), we evidenced an asymmetry between prograde and retrograde-planar components down to the most metal-poor values. We confirmed that stars on prograde-planar orbits are present well below –1.7 dex, that is, 1.0 dex below the typical limit suggested for that kind of population (e.g. Bensby et al. 2014; Fuhrmann et al. 2017). This prograde-planar subset represents 10% of the PGS VMP population, in the rotated action space of the decontaminated PGS between 0.2 and 0.7 (see Figs. 11 and E1). Comparably to Sestito et al. (2020), we observed traces of the VMP highly prograde (V<sub>φ</sub> > 180 km.s<sup>-1</sup>) and planar (Z<sub>max</sub> < 1.5 kpc) sample (corresponding to a subset of 2% of the total VMP sample) down to –2.9 dex with a 2σ confidence interval, in a sample four times larger. Below a metallicity of –3 dex, our analysis of the PGS VMP sample is limited simultaneously by the lack of reference spectroscopic values in APOGEE, LAMOST, GALAH, and GSP-spec needed to properly perform the decontamination, less precise PGS metallicities, and overall low statistics to detect strong signal.

The main findings from the orbital characteristics of this population can be summarised as follows:

- the GUMS thin disc in the action space (Fig. 8) and the Z<sub>max</sub> – R<sub>max</sub> distribution (Fig. 13) is essentially concentrated towards the Galactic plane and displays circular orbits, which greatly differs with the representation of our VMP highly prograde and planar population;

- the PGS VMP highly prograde and planar sample is more concentrated towards the inner MW and at larger Z than the higher-metallicity prograde-planar stars in the spatial

distribution (Fig. 12), compatible either with a short-scaled thick disc, or with the metal-poor inner Galaxy population of Rix et al. (2022) and Ardern-Arentsen et al. (2024);

- the PGS VMP highly prograde and planar sample is more concentrated towards the Galactic plane (Z<sub>max</sub> < 2 kpc) than the full VMP sample (33% versus 22%), while the classical GUMS halo shows significantly lower concentrations (14–18%). This calls for a contribution from a thick disc-like population, which in GUMS shows a similar concentration below Z<sub>max</sub> = 2 kpc (Fig. 13);

- the PGS VMP highly prograde and planar sample cannot solely be part of a canonical halo, as shown in the V<sub>φ</sub> distribution (residuals in the sixth panel of Fig. 14 ~ 150 km.s<sup>-1</sup>) and the action space (residuals in the bottom panels of Fig. 15 and Fig. E3). To remain agnostic on distribution functional forms and models, we assumed that the halo is a single component, with symmetrical prograde and retrograde tails of velocity distribution. When folding the retrograde over the prograde distribution of our sample at the lowest metallicities, we found that a stationary halo was not enough to suppress highly prograde VMP stars, while a prograde halo centred at 30 or 40 km.s<sup>-1</sup> potentially could.

These characteristics show that this population is a combination of some (possibly prograde) halo contribution and some thick disc-like component, reminiscent of the metal-weak thick disc (Morrison et al. 1990; Kordopatis et al. 2013, 2017; Beers et al. 2014; Carollo et al. 2019), or the Atari disc of Mardini et al. (2022). Different origins can be considered for this population. Using ARTEMIS simulations of MW-mass galaxies, Dillamore et al. (2024) examine the formation of the MW disc. In fact, they show that an early spin-up anti-correlates with the fraction of accreted stars and correlates with high halo mass. Since the MW has an early spin-up, most stars that formed the disc were most likely already inside the MW at the epoch of spin-up. Within the scope of an early spin-up, Semenov et al. (2018) and show that the 10% of their MW analogues that formed discs at an early epoch (more than 10 Gyr ago) had a lower metallicity limit located between –1 and –1.5 dex. We note, however, that Sotillo-Ramos et al. (2023) use the same suite of simulations to constrain the orbital properties of MW-like galaxies, solely based on their circularity, their stellar mass and their environment, and find that 20% of stars with [Fe/H] < –2.0 belong to kinematically cold geometrical thin discs, with some UMP analogues reaching ~50%. They also find that such stars were formed ex situ. In our case, the stars we consider to be disc-like (highly prograde and planar) make up for 2% of the PGS VMP sample. Such simulations, in light of our results, suggest that we might be just a few instruments away from uncovering a larger population of these singular stars. Belokurov & Kravtsov (2022) work with APOGEE and Gaia EDR3 data and evidenced Aurora, a low-metallicity in situ component with [Fe/H] ≤ –1.3 and mild net rotation; such population is in agreement with the prograde halo of Zhang et al. (2023), who made use of the XGBoost metallicity sample of Andrae et al. (2023). The test-particle simulations of Li et al. (2024a), created with Gaia RVS and Pristine DR1 data, reveal that a rapidly decelerating bar can transfer inner Galaxy stars towards very prograde-planar orbits, giving an explanation for the presence of stars with very low metallicities along the Galactic disc. With the same method, this time only using Gaia RVS data, Dillamore et al. (2023) find ridges in phase space interpreted as signatures of bar resonances, causing the formation of substructure in the disc. Ardern-Arentsen et al. (2024), using spectroscopic data from the
Pristine Inner Galaxy Survey (PIGS, Arentsen et al. 2020) and distances derived with StarHorse (Queiroz et al. 2018), suggest that within $1.5 < R < 3.5$ kpc and below $[Fe/H] = -2.0$, there could be a mix of a stationary halo component (40%) and a prograde halo component (60%). Finally, Li et al. (2024b) study the N-body simulations of barred disc galaxies, with spinning and stationary dark matter halos, obtained with the hydrodynamical code GIZMO (Hopkins 2015). They propose a mechanism to drive metal-poor stars onto disc orbits consisting of cooling the vertical oscillations of halo stars, leading to a flattening of $Z_{max}$ and $R_{max}$. Another possibility could involve accretion, as discussed in Re Fiorentin et al. (2021). One of the substructures identified, dubbed Icarus, is quite similar to the population investigated in this paper; nonetheless, it could be the remnant of a dwarf galaxy progenitor located on a prograde and low-inclination orbit, according to the study of its chemodynamics, in particular α-elements, metallicity and eccentricity. Simulations such as that of Sestito et al. (2021) point out the possibility of later accretions of disrupting galaxies, permeating the disc with low-metallicity stars.

To conclude, we must derive elemental abundances to disentangle the possible origins of the PGS VMP prograde-planar population. Hawkins et al. (2015) and Feltzing & Feuillet (2023), among other studies, show how powerful chemical tagging is, in particular when investigating the $[\alpha/Fe]-[Fe/H]$ and the Mg-Mn-Al-Fe planes, on the one hand to locate the canonical Galactic components, and on the other hand to effectively separate in situ from accreted stars. For example, the work of Fernández-Alvar et al. (2024) (see their Figs. 4 and 5) constrains the thin disc at $-1.0 < [Fe/H] < -0.7$ and shows that high $V_\phi$ stars with $[Fe/H] < -1.0$ have trends typical of early chemical enrichment, in opposition with thick disc trends. Other studies such as that of Dovgal et al. (2024) underline the importance of clarifying the chemical enrichment history of metal-poor disc stars, paving the way towards efficient use of upcoming spectroscopic facilities such as WEAVE (Jin et al. 2023), 4MOST (de Jong et al. 2019), DESI (DESI Collaboration 2016), and MOONS (Cirasuolo et al. 2020). Future prospects involve the spectral analysis of a statistically robust sample composed of canonical thin and thick disc stars as well as Pristine metal-poor and VMP stars, to confirm their metallicities, analyse, and compare their chemical content. We believe that this may be the missing piece needed to reconstruct a clear history of the Galactic disc.

**Acknowledgements.** We thank the anonymous referee for their insightful comments that helped clarify the contents of this paper. We warmly thank Pedro A. Palicio for useful discussions that helped improve the quality of this paper. NFM, VH, GK and IGR gratefully acknowledge support from the French National Research Agency (ANR) funded project "Pristine" (ANR-18-CE31-0017). GK, VH, FG and IGR gratefully acknowledge support from the French National Research Agency (ANR) funded project MWDisc (ANR-20-CE31-0004). AAA acknowledges support from the Herchel Smith Fellowship at the University of Cambridge and a Fitzwilliam College research fellowship supported by the Isaac Newton Trust. GT and GB acknowledge support from the Agencia Estatal de Investigación del Ministerio de Ciencia en Innovación (AEI-MICIN) and the European Regional Development Fund (ERDF) under grant number PID2020-118778GB-100/10.13039/501100011033 and the AEI under grant number CEX2019-000920-S. ES acknowledges funding through VIDI grant "Pushing Galactic Archaeology to its limits" (with project number VI.Vidi.193.093) which is funded by the Dutch Research Council (NWO). This research has been partially funded from a Spinoza award by NWO (SPI 78-411). This research was supported by the International Space Science Institute (ISSI) in Bern, through ISSI International Team project 540 (The Early Milky Way). SV thanks ANID (Beca Doctorado Nacional, folio 21220489) and the Millennium Nucleus ERIS (ERIS NCN2021017). This work has made use of data from the European Space Agency (ESA) mission Gaia

(https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. Funding for the Sloan Digital Sky Survey IV has been provided by the Alfred P. Sloan Foundation, the U.S. Department of Energy Office of Science, and the Participating Institutions. SDSS-IV acknowledges support and resources from the Center for High Performance Computing at the University of Utah. The SDSS website is www.sdss4.org. SDSS-IV is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS Collaboration including the Brazilian Participation Group, the Carnegie Institution for Science, Carnegie Mellon University, Center for Astrophysics | Harvard & Smithsonian, the Chilean Participation Group, the French Participation Group, Instituto de Astrofísica de Canarias, The Johns Hopkins University, Kavli Institute for the Physics and Mathematics of the Universe (IPMU) / University of Tokyo, the Korean Participation Group, Lawrence Berkeley National Laboratory, Leibniz Institut für Astrophysik Potsdam (AIP), Max-Planck-Institut für Astronomie (MPIA Heidelberg), Max-Planck-Institut für Astrophysik (MPA Garching), Max-Planck-Institut für Extraterrestrische Physik (MPE), National Astronomical Observatories of China, New Mexico State University, New York University, University of Notre Dame, Observatório Nacional / MCTI, The Ohio State University, Pennsylvania State University, Shanghai Astronomical Observatory, United Kingdom Participation Group, Universidad Nacional Autónoma de México, University of Arizona, University of Colorado Boulder, University of Oxford, University of Portsmouth, University of Utah, University of Virginia, University of Washington, University of Wisconsin, Vanderbilt University, and Yale University. This work made use of the Third Data Release of the GALAH Survey (Buder et al. 2021). The GALAH Survey is based on data acquired through the Australian Astronomical Observatory, under programs: A/2013B/13 (The GALAH pilot survey); A/2014A/25, A/2015A/19, A2017A/18 (The GALAH survey phase 1); A2018A/18 (Open clusters with HERMES); A2019A/1 (Hierarchical star formation in Ori OB1); A2019A/15 (The GALAH survey phase 2); A/2015B/19, A/2016A/22, A/2016B/10, A/2017B/16, A/2018B/15 (The HERMES-TESS program); and A/2015A/3, A/2015B/1, A/2015B/19, A/2016A/22, A/2016B/12, A/2017A/14 (The HERMES K2-follow-up program). We acknowledge the traditional owners of the land on which the AAT stands, the Gamilaray people, and pay our respects to elders past and present. This paper includes data that has been provided by AAO Data Central (datacentral.org.au). Guoshoujing Telescope (the Large Sky Area Multi-Object Fiber Spectroscopic Telescope LAMOST) is a National Major Scientific Project built by the Chinese Academy of Sciences. Funding for the project has been provided by the National Development and Reform Commission. LAMOST is operated and managed by the National Astronomical Observatories, Chinese Academy of Sciences.

## References

Abazajian, K., Adelman-McCarthy, J. K., Agüeros, M. A., et al. 2003, AJ, 126, 2081

Abdurro'uf, Accetta, K., Aerts, C., et al. 2022, ApJS, 259, 35

Aguado, D. S., Youakim, K., González Hernández, J. I., et al. 2019, MNRAS, 490, 2241

Andrae, R., Rix, H.-W., & Chandra, V. 2023, ApJS, 267, 8

Ardern-Arentsen, A., Monari, G., Queiroz, A. B. A., et al. 2024, MNRAS, 530, 3391

Arentsen, A., Starkenburg, E., Martin, N. F., et al. 2020, MNRAS, 491, L11

Bailer-Jones, C. A. L., Rybizki, J., Fouesneau, M., Demleitner, M., & Andrae, R. 2021, AJ, 161, 147

Beers, T. C., & Christiebly, N. 2005, ARA&A, 43, 531

Beers, T. C., Preston, G. W., & Shectman, S. A. 1985, AJ, 90, 2089

Beers, T. C., Preston, G. W., & Shectman, S. A. 1992, AJ, 103, 1987

Beers, T. C., Norris, J. E., Placco, V. M., et al. 2014, ApJ, 794, 58

Bellazzini, M., Massari, D., Ceccarelli, E., et al. 2024, A&A, 683, A136

Belokurov, V., & Kravtsov, A. 2022, MNRAS, 514, 689

Belokurov, V., Sanders, J. L., Fattahi, A., et al. 2020, MNRAS, 494, 3880

Bennett, M., & Bovy, J. 2019, MNRAS, 482, 1417

Bensby, T., Alves-Brito, A., Oey, M. S., Yong, D., & Meléndez, J. 2011, ApJ, 735, L46

Bensby, T., Feltzing, S., & Oey, M. S. 2014, A&A, 562, A71

Binney, J. 2012a, MNRAS, 426, 1324

Binney, J. 2012b, MNRAS, 426, 1328

Bond, H. E. 1970, ApJS, 22, 117

Bond, H. E. 1980, ApJS, 44, 517

Bonifacio, P., Caffau, E., Sestito, F., et al. 2019, MNRAS, 487, 3797

Bovy, J. 2015, ApJS, 216, 29

Bovy, J., & Rix, H.-W. 2013, ApJ, 779, 115

Bovy, J., Rix, H.-W., & Hogg, D. W. 2012, ApJ, 751, 131
Bressan, A., Marigo, P., Girardi, L., et al. 2012, MNRAS, 427, 127

Buder, S., Lind, K., Ness, M. K., et al. 2019, A&A, 624, A19

Buder, S., Sharma, S., Kos, J., et al. 2021, MNRAS, 506, 150

Caffau, E., Lombardo, L., Mashonkina, L., et al. 2023, MNRAS, 518, 3796

Carollo, D., Chiba, M., Ishigaki, M., et al. 2019, ApJ, 887, 22

Carollo, D., Christlieb, N., Tissera, P. B., & Sillero, E. 2023, ApJ, 946, 99

Chen, Y. H., Li, G. W., & Shu, H. 2022, Res. Astron. Astrophys., 22, 055008

Chiba, M., & Beers, T. C. 2000, AJ, 119, 2843

Chiba, M., & Yoshii, Y. 1998, AJ, 115, 168

Christlieb, N., Gustafsson, B., Korn, A. J., et al. 2004, ApJ, 603, 708

Cirasuolo, M., Fairley, A., Rees, P., et al. 2020, Messenger, 180, 10

Comparat, J., Delubac, T., Jouvel, S., et al. 2016, A&A, 592, A121

Cordoni, G., Da Costa, G. S., Yong, D., et al. 2021, MNRAS, 503, 2539

Cropper, M., Katz, D., Sartoretti, P., et al. 2018, A&A, 616, A5

Cui, X.-Q., Zhao, Y.-H., Chu, Y.-Q., et al. 2012, Res. Astron. Astrophys., 12, 1197

Dawson, K. S., Schlegel, D. J., Ahn, C. P., et al. 2013, AJ, 145, 10

Deason, A. J., Belokurov, V., & Evans, N. W. 2011, MNRAS, 411, 1480

Deason, A. J., Belokurov, V., Koposov, S. E., et al. 2017, MNRAS, 470, 1259

de Jong, R. S., Agertz, O., Berbel, A. A., et al. 2019, Messenger, 175, 3

De Silva, G. M., Freeman, K. C., Bland-Hawthorn, J., et al. 2015, MNRAS, 449, 2604

DESI Collaboration (Aghamousa, A., et al.) 2016, arXiv e-prints [arXiv:1611.00036]

Dillamore, A. M., Belokurov, V., Evans, N. W., & Davies, E. Y. 2023, MNRAS, 524, 3596

Dillamore, A. M., Belokurov, V., Kravtsov, A., & Font, A. S. 2024, MNRAS, 527, 7070

Di Matteo, P., Spite, M., Haywood, M., et al. 2020, A&A, 636, A115

Dovgal, A., Venn, K. A., Sestito, F., et al. 2024, MNRAS, 527, 7810

Eggen, O. J., Lynden-Bell, D., & Sandage, A. R. 1962, ApJ, 136, 748

Eilers, A.-C., Hogg, D. W., Rix, H.-W., & Ness, M. K. 2019, ApJ, 871, 120

El-Badry, K., Bland-Hawthorn, J., Wetzel, A., et al. 2018, MNRAS, 480, 652

Feltzing, S., & Feuillet, D. 2023, ApJ, 953, 143

Fernández-Alvar, E., Kordopatis, G., Hill, V., et al. 2021, MNRAS, 508, 1509

Fernández-Alvar, E., Kordopatis, G., Hill, V., et al. 2024, A&A, 685, A151

Re Fiorentin, P., Spagna, A., Lattanzi, M. G., & Cignoni, M. 2021, ApJ, 907, L16

Frebel, A., & Norris, J. E. 2013, in Planets, Stars and Stellar Systems, eds. T. D. Oswalt, & G. Gilmore, Galactic Structure and Stellar Populations, 5, 55

Frebel, A., & Norris, J. E. 2015, ARA&A, 53, 631

Freeman, K., & Bland-Hawthorn, J. 2002, ARA&A, 40, 487

Fuhrmann, K., Chini, R., Kaderhandt, L., & Chen, Z. 2017, MNRAS, 464, 2610

Gaia Collaboration (Brown, A. G. A., et al.) 2016a, A&A, 595, A2

Gaia Collaboration (Prusti, T., et al.) 2016b, A&A, 595, A1

Gaia Collaboration (Vallenari, A., et al.) 2023a, A&A, 674, A1

Gaia Collaboration (Recio-Blanco, A., et al.) 2023b, A&A, 674, A38

Gaia Collaboration (Montegriffo, P., et al.) 2023c, A&A, 674, A33

Gallart, C., Surot, F., Cassisi, S., et al. 2024, A&A, 687, A168

García Pérez, A. E., Allende Prieto, C., Holtzman, J. A., et al. 2016, AJ, 151, 144

Gilmore, G., & Reid, N. 1983, MNRAS, 202, 1025

GRAVITY Collaboration (Abuter, R., et al.) 2020, A&A, 636, L5

Hawkins, K., Jofré, P., Masseron, T., & Gilmore, G. 2015, MNRAS, 453, 758

Hayden, M. R., Holtzman, J. A., Bovy, J., et al. 2014, AJ, 147, 116

Hayden, M. R., Bovy, J., Holtzman, J. A., et al. 2015, ApJ, 808, 132

Haywood, M., Di Matteo, P., Lehnert, M. D., et al. 2018, ApJ, 863, 113

Hegedüs, V., Mészáros, S., Jofré, P., et al. 2023, A&A, 670, A107

Hopkins, P. F. 2015, MNRAS, 450, 53

Hunt, E. L., & Reffert, S. 2023, A&A, 673, A114

Ibata, R. A., Gilmore, G., & Irwin, M. J. 1994, Nature, 370, 194

Jin, S., Trager, S. C., Dalton, G. B., et al. 2023, MNRAS, 530, 2688

Katz, D., Sartoretti, P., Guerrier, A., et al. 2023, A&A, 674, A5

Keller, S. C., Schmidt, B. P., Bessell, M. S., et al. 2007, PASA, 24, 1

Kielty, C. L., Venn, K. A., Sestito, F., et al. 2021, MNRAS, 506, 1438

Koppelman, H. H., Hagen, J. H. J., & Helmi, A. 2021, A&A, 647, A37

Kordopatis, G., Recio-Blanco, A., de Laverny, P., et al. 2011, A&A, 535, A107

Kordopatis, G., Gilmore, G., Wyse, R. F. G., et al. 2013, MNRAS, 436, 3231

Kordopatis, G., Wyse, R. F. G., Gilmore, G., et al. 2015, A&A, 582, A122

Kordopatis, G., Wyse, R. F. G., Chiappini, C., et al. 2017, MNRAS, 467, 469

Kordopatis, G., Schultheis, M., McMillan, P. J., et al. 2023, A&A, 669, A104

Lardo, C., Mashonkina, L., Jablonka, P., et al. 2021, MNRAS, 508, 3068

Lee, Y. S., Beers, T. C., Sivarani, T., et al. 2008, AJ, 136, 2022

Li, H., Tan, K., & Zhao, G. 2018, ApJS, 238, 16

Li, C., Yuan, Z., Monari, G., et al. 2024a, A&A, 690, A26

Li, X., Shlosman, I., Pfenniger, D., & Heller, C. 2024b, ApJ, 971, L36

Lindegren, L., Klioner, S. A., Hernández, J., et al. 2021a, A&A, 649, A2

Lindegren, L., Bastian, U., Biermann, M., et al. 2021b, A&A, 649, A4

Lucchesi, R., Lardo, C., Jablonka, P., et al. 2022, MNRAS, 511, 1004

Lucey, M., Hawkins, K., Ness, M., et al. 2021, MNRAS, 501, 5981

Luo, A. L., Zhao, Y.-H., Zhao, G., et al. 2015, Res. Astron. Astrophys., 15, 1095

Mackereth, J. T., & Bovy, J. 2018, PASP, 130, 114501

Mackereth, J. T., Bovy, J., Leung, H. W., et al. 2019, MNRAS, 489, 176

Majewski, S. R., Schiavon, R. P., Frinchaboy, P. M., et al. 2017, AJ, 154, 94

Mardini, M. K., Frebel, A., Chiti, A., et al. 2022, ApJ, 936, 78

Martin, N. F., Starkenburg, E., Yuan, Z., et al. 2024, A&A, 692, A115

Matsuno, T., Starkenburg, E., Balbinot, E., & Helmi, A. 2024, A&A, 685, A59

McMillan, P. J. 2017, MNRAS, 465, 76

Mendes de Oliveira, C., Ribeiro, T., Schoenell, W., et al. 2019, MNRAS, 489, 241

Morrison, H. L., Flynn, C., & Freeman, K. C. 1990, AJ, 100, 1191

Nepal, S., Chiappini, C., Queiroz, A. B. A., et al. 2024, A&A, 688, A167

Norris, J., Bessell, M. S., & Pickles, A. J. 1985, ApJS, 58, 463

Perryman, M. A. C., de Boer, K. S., Gilmore, G., et al. 2001, A&A, 369, 339

Pillepich, A., Springel, V., Nelson, D., et al. 2018, MNRAS, 473, 4077

Pöder, S., Benito, M., Pata, J., et al. 2023, A&A, 676, A134

Queiroz, A. B. A., Anders, F., Santiago, B. X., et al. 2018, MNRAS, 476, 2556

Randich, S., Gilmore, G., Magrini, L., et al. 2022, A&A, 666, A121

Recio-Blanco, A., de Laverny, P., Allende Prieto, C., et al. 2016, A&A, 585, A93

Recio-Blanco, A., de Laverny, P., Palicio, P. A., et al. 2023, A&A, 674, A29

Reid, M. J., & Brunthaler, A. 2020, ApJ, 892, 39

Rix, H.-W., Chandra, V., Andrae, R., et al. 2022, ApJ, 941, 45

Robin, A. C., Reylé, C., Derrière, S., & Picaud, S. 2003, A&A, 409, 523

Robin, A. C., Luri, X., Reylé, C., et al. 2012, A&A, 543, A100

Robin, A. C., Bienaymé, O., Fernández-Trincado, J. G., & Reylé, C. 2017, A&A, 605, A1

Roman, N. G. 1950, ApJ, 112, 554

Ruchti, G. R., Fulbright, J. P., Wyse, R. F. G., et al. 2011, ApJ, 737, 9

Santistevan, I. B., Wetzel, A., Sanderson, R. E., et al. 2021, MNRAS, 505, 921

Schlegel, D. J., Finkbeiner, D. P., & Davis, M. 1998, ApJ, 500, 525

Schuster, W. J., & Nissen, P. E. 1989, A&A, 222, 69

Schwarzschild, M., & Schwarzschild, B. 1950, ApJ, 112, 248

Searle, L., & Zinn, R. 1978, ApJ, 225, 357

Semenov, V. A., Conroy, C., Chandra, V., Hernquist, L., & Nelson, D. 2024, ApJ, 962, 84

Sestito, F., Longeard, N., Martin, N. F., et al. 2019, MNRAS, 484, 2166

Sestito, F., Martin, N. F., Starkenburg, E., et al. 2020, MNRAS, 497, L7

Sestito, F., Buck, T., Starkenburg, E., et al. 2021, MNRAS, 500, 3750

Sestito, F., Hayes, C. R., Venn, K. A., et al. 2024, MNRAS, 528, 4838

Sofue, Y. 2015, PASJ, 67, 75

Sotillo-Ramos, D., Bergemann, M., Friske, J. K. S., & Pillepich, A. 2023, MNRAS, 525, L105

Starkenburg, E., Oman, K. A., Navarro, J. F., et al. 2017a, MNRAS, 465, 2122

Starkenburg, E., Martin, N., Youakim, K., et al. 2017b, MNRAS, 471, 2587

Starkenburg, E., Aguado, D. S., Bonifacio, P., et al. 2018, MNRAS, 481, 3838

Steinmetz, M., Matijevič, G., Enke, H., et al. 2020, AJ, 160, 82

Thomas, G. F., Laporte, C. F. P., McConnachie, A. W., et al. 2019, MNRAS, 483, 3119

Tian, H., Liu, C., Xu, Y., & Xue, X. 2019, ApJ, 871, 184

Tsantaki, M., Pancino, E., Marrese, P., et al. 2022, A&A, 659, A95

Tumlinson, J. 2010, ApJ, 708, 1398

Vasiliev, E., & Baumgardt, H. 2021, MNRAS, 505, 5978

Venn, K. A., Kielty, C. L., Sestito, F., et al. 2020, MNRAS, 492, 3241

Vieira, K., Carraro, G., Korchagin, V., et al. 2022, ApJ, 932, 28

Viswanathan, A., Starkenburg, E., Matsuno, T., et al. 2024a, A&A, 683, L11

Viswanathan, A., Yuan, Z., Ardern-Arentsen, A., et al. 2024b, A&A, submitted [arXiv:2405.13124]

White, S. D. M., & Springel, V. 2000, in The First Stars, eds. A. Weiss, T. G. Abel, & V. Hill, 327

Wisotzki, L., Koehler, T., Groote, D., & Reimers, D. 1996, A&AS, 115, 227

Wisotzki, L., Christlieb, N., Bade, N., et al. 2000, A&A, 358, 77

Xiang, M., & Rix, H.-W. 2022, Nature, 603, 599

Yuan, Z., Li, C., Martin, N. F., et al. 2024, A&A, 691, L1

Zhang, H., Ardern-Arentsen, A., & Belokurov, V. 2023, MNRAS, 533, 889

Zolotov, A., Willman, B., Brooks, A. M., et al. 2009, ApJ, 702, 1058
# Appendix A: PGS quality selection path

Table A.1 describes the chronological cuts applied to PGS in order to build the sample used for this work. As a reminder, the published PGS contains ~ 52 million stars; our selection represents ~ 6 % of the original catalogue, corresponding to roughly 3 million stars.

# Appendix B: Isochrone filtering: GALAH, LAMOST, and GSP-spec

In the same fashion as in Sec. 4.1, we briefly discuss the validation of the isochrone filtering method for GALAH, LAMOST and GSP-spec. Table B1 summarises the sample size of each cross-match between PGS and spectroscopic surveys before and after applying the isochrone filtering method.

GALAH is the sample for which we have the lowest statistics. However, we are able to evidence in the left panel of Fig. B.1a the same horizontal sequence of outliers identified in the PGS-APOGEE cross-correlation. This sequence makes up for 13.5 % of the total number of PGS-GALAH metal-poor candidates. After filtering, 63 % of GALAH interlopers are removed. Although this score is less satisfying than with APOGEE, most of the horizontal sequence is effectively removed as seen in the top right panel of Fig. B.1a. Additionally, 5.8 % of PGS-GALAH metal-poor sample are also discarded.

The LAMOST sample exhibits a bold horizontal sequence of interlopers in the left panel of Fig. B.2a. As seen in the left panel of Fig. B.2b, they are for the most part dwarf stars, which are removed as part of our selection function; nonetheless, the clump at $M_G \sim 0$ similar to the PGS-APOGEE cross-correlation is still present, but can be effectively removed with the isochrone filtering method. The interlopers make up for 11.2 % of the PGS-LAMOST metal-poor subsample. After filtering, 89 % of LAMOST interlopers are discarded, with an additional 56.4 % of the true PGS-LAMOST metal-poor stars (mostly dwarf stars) and effective removal of the interloper sequence in the right panel of Fig. B.2a.

The GSP-spec sample, in addition to the horizontal interloper sequence, displays a vertical interloper sequence towards solar [Fe/H] in the left panel of Fig. B3a, this time corresponding to metal-poor (most likely spectroscopic) interlopers. Interlopers are clumped in the reddest part of the CMD and are mainly giant stars, see left panel of Fig. B3b. They make up for 31 % of the PGS - GSP-spec metal-poor subsample. The isochrone filtering method discards both interloper sequences quite successfully. After filtering, 67 % of GSP-spec interlopers are discarded, with an additional 2.2 % of the true PGS - GSP-spec metal-poor stars.

<figure><img src="image_18.png" /><figcaption>**Fig. B.1.** PGS - GALAH before (left) and after (right) applying the isochrone filtering method. Top: $[Fe/H]_{GALAH}$ versus $[Fe/H]_{PGS}$ (Fig. B.1a). Bottom: CMD (Fig. B.1b).</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>**Fig. B.2.** PGS - LAMOST before (left) and after (right) applying the isochrone filtering method. Top: $[Fe/H]_{LAMOST}$ versus $[Fe/H]_{PGS}$ (Fig. B.2a). Bottom: CMD (Fig. B.2b).</figcaption></figure>

# Appendix C: GSP-spec flags quality selection

Table C1 describes the flag and stellar parameter selection of our GSP-spec sample.

# Appendix D: On the choice of spectroscopic contaminants

To decide which surveys are the best fit for the spectroscopic decontamination of PGS, we investigated in Fig. D2 the $V_\phi$ distribution in the same [Fe/H] intervals defined in Fig. 6, and compared APOGEE, GALAH, LAMOST and GSP-spec. We refer the reader to Sec. 4.2 for a complete description of the procedure used to isolate spectroscopic contaminants, and Fig. D1 to

visualise the effective selection relative to GALAH, LAMOST and GSP-spec.

Between -1.0 and -1.7, and between -1.7 and -4.0, GSP-spec contaminants follow a trend distinct from the other spectroscopic surveys. In fact, below -1.7, the bulk of the contamination lies at ~ 250 km.s⁻¹, whereas the contamination of the other surveys is centred around 0 km.s⁻¹ and spreads down to retrograde regions. In particular, GSP-spec contamination exceeds true PGS counts, which drives us to identify the cause for this overestimation.

We suspect that GSP-spec is incomplete in the metal-poor domain, due to the data selection itself; strictly following the quality flags advised in Table 2 of Recio-Blanco et al. (2023) can remove a significant number of metal-poor stars, due to
<table><caption>**Table A.1.** Description of the various cuts applied to the *Pristine-Gaia* synthetic catalogue, resulting in the version used for this work.</caption><thead><tr><th>Cut</th><th>Number of stars</th><th>Percentage removed from PGS (%)</th></tr></thead><tbody><tr><td>No cut (published PGS)</td><td>52 300 000</td><td>0</td></tr><tr><td>Gaia cuts except gaiadr3.vari_summary</td><td>41 612 344</td><td>20.4</td></tr><tr><td>(BP - RP)<sub>0</sub> > 0.5 mag and mcfrac_CaHKsyn > 0 (no MC draw available)</td><td>39 446 013</td><td>4.1</td></tr><tr><td>M<sub>G</sub> < 7.5</td><td>39 217 953</td><td>0.44</td></tr><tr><td>mcfrac_CaHKsyn > 0.9</td><td>37 069 113</td><td>4.11</td></tr><tr><td>Sources with available RVS orbits</td><td>12 777 188</td><td>46.45</td></tr><tr><td>[Fe/H] > -4.0</td><td>12 776 653</td><td>0.001</td></tr><tr><td>E(B - V) < 0.5</td><td>12 737 484</td><td>0.08</td></tr><tr><td>Pvar < 0.3 and d_CaHK < 0.05</td><td>9 548 136</td><td>6.10</td></tr><tr><td>gaiadr3.vari_summary</td><td>9 450 368</td><td>0.19</td></tr><tr><td>Catalogue of globular and open cluster stars</td><td>9 432 286</td><td>0.03</td></tr><tr><td>Dwarf and giant separation + isochrone filtering</td><td>2 880 338</td><td>12.53</td></tr></tbody></table>

<table><caption>**Table B1.** Sample size of the cross-matches between PGS and relevant spectroscopic catalogues.</caption><thead><tr><th>Survey</th><th>Before filtering (giants + dwarfs)</th><th>After filtering (only giants)</th></tr></thead><tbody><tr><td>APOGEE</td><td>178 637</td><td>92 765</td></tr><tr><td>GALAH</td><td>42 424</td><td>14 424</td></tr><tr><td>LAMOST</td><td>2 002 355</td><td>313 833</td></tr><tr><td>Gaia RVS (GSP-spec)</td><td>2 078 139</td><td>1 101 704</td></tr></tbody></table>

**Notes.** The cross-matches are used for the validation of the isochrone filtering method described in Sec. 4.1.

<figure><img src="image_20.png" /><figcaption>**Fig. B3.** PGS - GSP-spec before (left) and after (right) applying the isochrone filtering method. Top: $[Fe/H]_{GSP-spec}$ versus $[Fe/H]_{PGS}$ (Fig. B3a). Bottom: CMD (Fig. B3b).</figcaption></figure>

lower accuracy in the retrieval of their atmospheric parameters and $[Fe/H]$ (quality of the spectra). However, to be able to populate the contamination throughout the entire $V_\phi$ space, the cuts applied to our GSP-spec sample were loosened (see Table C1). Assuming that PGS is complete for halo (retrograde)

<table><caption>**Table C1.** Flag and stellar parameter selection of our GSP-spec sample.</caption><tr><td>vbroadT</td><td>0,1</td></tr><tr><td>vradT</td><td>0,1</td></tr><tr><td>vbroadG</td><td>0,1</td></tr><tr><td>vradG</td><td>0,1</td></tr><tr><td>vbroadM</td><td>0,1</td></tr><tr><td>vradM</td><td>0,1</td></tr><tr><td>negFlux</td><td>0</td></tr><tr><td>KMgiantPar</td><td>0</td></tr><tr><td>extrapol</td><td>0,1,2</td></tr><tr><td>fluxNoise</td><td>0,1,2</td></tr><tr><td>nullFluxErr</td><td>0,1</td></tr><tr><td>nanFlux</td><td>0,1</td></tr><tr><td>emission</td><td>0,1</td></tr><tr><td>T<sub>eff</sub> (K)</td><td>[3800, 6000]</td></tr><tr><td>RUWE</td><td>[0, 1.4]</td></tr><tr><td>log g (dex)</td><td>[0, 5]</td></tr></table>

**Notes.** For a more complete description of the flags, see Table 2 of Recio-Blanco et al. (2023).

bright stars (for which the Pristine pipeline performs best), we checked the completeness of APOGEE, GALAH, LAMOST and GSP-spec with respect to PGS, as a function of PGS metallicity in Fig. D3. We selected sources with $V_\phi < 80$ km.s$^{-1}$ and $phot_g_mean_mag < 12$ mag (to include the bulk of GSP-spec stars, see Fig. 2 of Recio-Blanco et al. 2023).

We observe overall constant completeness with decreasing metallicity for APOGEE and LAMOST; still, the completeness of GSP-spec decreases with metallicity. Below -2.0 dex, it is downsized by a factor 3. Similarly to what is seen for retrograde stars, GSP-spec metal-poor sources might be missing in the high $V_\phi$ regions, leading to an overrepresentation of metal-rich contamination dominating the sample. Figure D4 confirms that simply rescaling GSP-spec contaminants by that factor enables to recover values coherent with other spectroscopic trends as well as PGS itself, and that GSP-spec is incomplete at very low metallicities. For that reason, we did not use GSP-spec to decontaminate PGS.

Regarding GALAH contaminants, it is complex to infer a trend for completeness because it is almost constantly null, due to very small statistics available after cross-matching with PGS. Hence, we did not use GALAH to decontaminate PGS.
<figure><img src="image_21.png" /><figcaption>Fig. D1. $\Delta_{[Fe/H]} = [Fe/H]_{PGS} - [Fe/H]_{survey}$ distribution (black). The Gaussian fit to the distribution is overplotted in red dashed lines. $\sigma$ and $2\sigma$ limits are respectively indicated in blue and teal vertical lines. The PGS-survey contaminants, i.e. sources with $\Delta_{[Fe/H]} = [Fe/H]_{PGS} - [Fe/H]_{survey} > 2\sigma$ are colour-coded in purple. Left panel: $-1.0 < [Fe/H]_{PGS} < -0.3$. Middle panel: $-1.7 < [Fe/H]_{PGS} < -1.0$. Right panel: $-4.0 < [Fe/H]_{PGS} < -1.7$. From top to bottom: comparison with respect to GALAH, LAMOST and GSP-spec.</figcaption></figure>

## Appendix E: Decontamination and chemodynamical analysis with LAMOST

We show in Fig. E1 the decontamination of PGS using LAMOST, similarly to the methodology of Sec. 5 with APOGEE. Figure E2 and Fig. E3 are relevant to the discussion of Sec. 6.2.
<figure><img src="image_22.png" /><figcaption>**Fig. D2.** $V_\phi$ distribution of PGS, together with spectroscopic contaminants scaled by the size of PGS in each [Fe/H] subsample. Top panel: $-1.0 < [\text{Fe/H}]_{\text{PGS}} < -0.3$. Middle panel: $-1.7 < [\text{Fe/H}]_{\text{PGS}} < -1.0$. Bottom panel: $-4.0 < [\text{Fe/H}]_{\text{PGS}} < -1.7$.</figcaption></figure>

<figure><img src="image_23.png" /><figcaption>**Fig. D3.** Completeness of halo stars for each spectroscopic survey as a function of [Fe/H].</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>**Fig. D4.** Same as Fig. D2, but only with GSP-spec. Magenta: Original GSP-spec contaminants. Green: GSP-spec contaminants downgraded by the decrease factor in completeness of halo stars (see Fig. D3). The error bar is divided by the square root of the decrease factor.</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>**Fig. E1.** Densities in the action space. The grid corresponds to areas in the parameter space delimited in ($J_z - J_r$) and $L_z$. From top to bottom: first panel: density of LAMOST contaminants. Second panel: PGS decontaminated pixel-by-pixel, using LAMOST. Third panel: residuals percentages of the pixel-by-pixel PGS decontamination using LAMOST. Fourth panel: PGS decontaminated by the mean of the LAMOST contamination.</figcaption></figure>
<figure><img src="image_26.png" /></figure>

**Fig. E2.** $V_\phi$ distribution for metallicities between -1.7 and -4.0. From top to bottom: panels 1 and 2: PGS decontamination by LAMOST. Panels 3-8: Folding of the retrograde counterparts over the prograde counterparts of the decontaminated PGS. The black dashed line indicates the location of our presumed halo. In panel 3, the folding is done assuming a theoretical halo. In the remaining panels, it is done assuming a prograde halo shifted by $10 \text{ km.s}^{-1}$ at every panel. The red dashed lines delimit the location of the population we aim to characterise. Same as Fig. 14 of Sec. 6.2.
<figure><img src="image_27.png" /></figure>

**Fig. E3.** Rotated action space for metallicities between -4.0 ≪ [Fe/H] < -1.7. The grey dashed lines correspond to different values of $L_z/J_{tot}$: [-0.95, -0.8, -0.5, 0, 0.5, 0.8, 0.95]. Left column: decontamination process and folding using LAMOST contamination. Top row: density of spectroscopic contaminants. Middle row: pixel-by-pixel decontamination. Bottom row: folding of the retrograde components over the prograde components, for PGS decontaminated pixel-by-pixel. In this space, the folding axis of symmetry is the grey dashed line at $L_z/J_{tot} = 0$. Right column: same as the left panels, but the three actions were computed shifting the observed $V_\phi$ by 40 km.s$^{-1}$. Same as Fig. 15 of Sec. 6.2.