Open Access
Issue
A&A
Volume 688, August 2024
Article Number A147
Number of page(s) 9
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202347450
Published online 15 August 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

Dark matter (DM) is central to the standard cosmological model (LCDM), providing gravitational support for the formation of galaxies and systems of galaxies (Mo et al. 2010). Its existence is backed by a plethora of observational data, including galaxies’ rotation curves (Zwicky 1933; Rubin et al. 1970), strong and weak gravitational lensing effects (Massey et al. 2010; Clowe et al. 2006), and even the presence of baryonic acoustic oscillation (Planck Collaboration VI. 2020) in the earlier gravitational wells observed in the cosmic microwave background (CMB).

Despite these successes, we still lack a precise detection either in the laboratory (Bernabei et al. 2022; Aprile et al. 2023; Barberio et al. 2023; Amaré et al. 2022) or by indirect astrophysical observation (Abdalla et al. 2022; Acharyya et al. 2023; Abe et al. 2023). Many theoretical candidates arise out of physics beyond the standard model (BSMP), like weakly interacting massive particles (WIMPS). These massive (mDM ∼ 100 GeV) particles could weakly interact with nucleons, and therefore their signals have been looked at by several laboratory and accelerator experiments. Also, their annihilation signals could be detectable through γ-ray emission by high-energy telescopes. Despite a massive experimental effort, DM remains a theoretical hypothesis, albeit one with impressive empirical support.

Other DM candidates could be more massive, such as primordial black holes (PBHs) (Carr & Kühnel 2020), which were recently constrained with a series of consistency tests. Nowadays, there is still room to be an essential contributor to the DM content but this is limited to small windows in mass (Villanueva-Domingo et al. 2021). Other candidates include massive ultralight particles (ULDMs) that could reach masses as low as 10−23 eV (Hui et al. 2017). Therefore, the possible DM mass range remains unconstrained today. Additionally, DM particles could interact with themselves, and have or not have spin, and other properties, including their mass, remain elusive.

Depending on the nature of the DM particle, there are relevant changes in the structure and number of DM halos and subhalos (Zavala & Frenk 2019). For example, some candidates, like warm dark matter (WDM), introduce a cut-off scale in the initial power spectrum of mass fluctuations (mDM ∼ 1 KeV), and others a scale during the non-linear evolution phase in which the DM particles self-interact (SIDM) (Tulin & Yu 2018). Both processes change the abundance of DM subhalos and the density profiles of the DM halos in comparison to the predictions of the cold dark matter (CDM) model.

Buschmann et al. (2018) show, using CDM simulations, that the gravitational pull of DM subhalos affects the distribution of stars in galactic halos, and could be used to discover dark subhalos (those without star formation in situ, a precise prediction of the CDM model) and also allow one to test the nature of the DM particles itself. This work simulated a passing DM subhalo’s perturbation to the phase space stellar distribution. Stars are pulled towards the subhalo as it passes, leaving a distinctive feature in the halo stars’ velocity and number density, known as a wake. This phenomenon was previously analytically described by Weinberg (1986) due to the gravitational friction that provokes the orbital decay of the satellite galaxies that inhabit the DM subhalo. There have been some efforts to quantify the magnitude of phase-space perturbations caused by the passage of DM subhalos using simulations and their possible detection (Bazarov et al. 2022).

Garavito-Camargo et al. (2019) used CDM simulations to quantify the impact of the LMC’s passage on the density and kinematics of the Milky Way (MW) DM halo and the observability of these structures in the MW’s stellar halo. Their results indicated a pronounced wake, which could be decomposed in a transient response and a collective one in both the DM and stellar halo distributions. Such an effect was observed for the first time in the MW halo stars by Conroy et al. (2021). The authors studied the effects induced by the Magellanic Clouds System (MCS) merging in selected samples of MW stars with precise Gaia Satellite measurements. This detection, and the increasing availability of stellar data from Gaia DR3 release, has paved the way for a more precise measurement of the effect. Measurements of the wake on the perturber systems of reference will allow pursuit further testing of the DM particle’s nature using this merger data, as was proposed recently by Foote et al. (2023) (in the context of systems of galaxies see Furlanetto & Loeb 2002; Buehler & Desjacques 2023). Furthermore, Aguilar-Argüello et al. (2022) and Cunningham et al. (2020) conducted studies on the decomposition into spherical harmonics of both density and velocity, respectively, in order to quantify the response of the DM halo to the passage of the Large Magellanic Cloud (LMC).

Based on the findings of Conroy et al. (2021), we have used Gaia Data Release 3 to study the DM subhalo of the Magellanic Clouds. Our code developments will also allow us to apply the methodology to other MW satellite galaxies and globular clusters and even to develop methods of detecting the presence of the dark subhalos predicted by the CDM model.

This work is organised as follows. In Sect. 2, we present the data samples selection methodology used to identify the effects of the DM subhalo of the LMC on the MW stellar halo. In Sect. 3, we briefly describe the theoretical model. Meanwhile, we present our results in Sect. 4. Section 5 presents the conclusions and future perspectives. Finally, we present in Appendix A a description of the coordinates’ transformation, in Appendix B the machine learning algorithm used to estimate radial velocities, and in Appendix C a validation of our method of inferring distances.

2. Data reduction

We studied the gravitational response of the MW’s halo to the passage of the LMC in its orbit. To achieve this, we used the data from Gaia Data Release 3 (Gaia Collaboration 2023, 2016) and created two catalogues of halo tracers; namely, K-giant stars and RR-Lyrae ones. We followed the steps proposed by Conroy et al. (2021) to address this task.

2.1. K-giant dataset

To construct the K-giant catalogue, we started the analysis with 162 240 774 sources characterised by ruwe values below 1.4, parallax measurements lower than 0.2 mas, and galactic latitudes of |b|> 10° to remove the galactic plane. To ensure data quality, we performed a series of cleaning procedures. First, we eliminated sources lacking proper motion and photometric data. To account for dust extinction, we obtained the dust map from Green (2018) and considered the Schlegel, Finkbeiner and Davis (SFD) map to derive the excess colour, E(B − V). Subsequently, we discarded all sources with E(B − V) > 0.3. To obtain the corrected magnitudes, we considered the following coefficients: AG/E(B − V) = 2.4, ABP/E(B − V) = 2.58, and ARP/E(B − V) = 1.65. To focus solely on the giant branch, we restricted the selection to sources satisfying the condition 1.4 < (BP* − RP*) < 2, where BP* and RP* represent the corrected magnitudes. Next, following Riello et al. (2021), we performed the 3σ cut upon the corrected BP and RP flux excess factor (C*). After completing the data-cleaning process, to ensure the purity of our catalogue specifically for K-giant stars, we performed a cross-match with the spectral types provided by Gaia (gaiadr3.astrophysical parameters) and obtained a dataset of 490 669 sources. Finally, we restricted our analysis to objects within a galactocentric distance between 30 kpc and 100 kpc, leaving 245 086 sources. Among them, only 10 989 had radial velocities measured by Gaia (Katz et al. 2023).

To estimate the radial velocity for the remaining sources, we employed a machine learning algorithm, specifically a RandomForestRegressor (Pedregosa et al. 2011). The accuracy of our model is 87.0% (see Appendix B for details).

To determine the photometric distance, we used the MIST code (Dotter 2016; Choi et al. 2016; Paxton et al. 2011) to generate an isochrone with the specific LMC parameters, which are an age of 10 Gyr and a metallicity of [Fe/H] = −1.5. We restricted the isochrone to an effective temperature from 3800 K to 4400 K, and fitted the polynomial equation

M G = 2.8894 ( B P R P ) 2 11.9263 ( B P R P ) + 8.7151 . $$ \begin{aligned} M_G = 2.8894(BP^{*}-RP^{*})^2-11.9263(BP^{*}-RP^{*})+8.7151. \end{aligned} $$(1)

To validate our distance inference method, we compared our calculated distances with reference values for some globular clusters, the Magellanic Clouds (LMC and SMC), and some satellite galaxies (see Appendix C for details). Our method successfully reproduced tabulated distances, since the mean calculated distances differed by less than 20%.

Afterwards, we implemented several masks to exclude known objects from our analysis. Specifically, we applied angular and distance masks to all the globular clusters listed in Harris (1996) and all the satellite galaxies reported in Drlica-Wagner et al. (2020).

Following the methodology proposed by Conroy et al. (2021), we employed proper motions to eliminate structures linked to the Sagittarius stream. To achieve this, we initially correct the proper motions due to the solar reflex motion, with the Gala package (Price-Whelan 2017; Price-Whelan et al. 2020). The used parameters were R = 8.122 kpc (Abuter 2019), (VR, ⊙, Vϕ, ⊙, VZ, ⊙) = (−12.9, 245.6, 7.78) km s−1 (Drimmel & Poggio 2018), and the distance of the Sun from the Galactic mid-plane, Z = 20.8 pc (Bennett & Bovy 2019). For b > 10° and |BSgr|< 15°, where BSgr is the latitude in the frame of the Sagittarius orbital plane, we removed part of the northern arm of the stream by taking out the sources that have μα* > −1.3 mas yr−1, −0.4 < μδ < 0.3 mas yr−1, and μδ > 1.7μα* + 0.4. To eliminate the rest of the north arm, we applied a mask to the region with coordinates b > 0° and 180° < l < 210°. The final selection of sources was based on proper motions; that is, we kept only those that satisfied μ α* 2 + ( μ δ +0.1) 2 < 0.5 2 $ \mu_{\alpha*}^2+(\mu_{\delta}+0.1)^2 < 0.5^2 $ (Conroy et al. 2021). By implementing this criterion, one effectively excludes disk stars, stars belonging to the Large and Small Magellanic Clouds, the Sagittarius dwarf spheroidal, and other Sagittarius arms. After this matching, our final dataset of K giants had 6058 sources.

2.2. RR-Lyrae dataset

To build the RR-Lyrae catalogue, we started the process by using the 271 779 sources catalogued as RR-Lyrae variables by Gaia. Initially, we performed data cleaning by keeping stars with ruwe < 1.4, and excluding those lacking metallicity, as well as those with errors in metallicity exceeding the absolute value of the metallicity itself. We then discarded all sources with E(B − V) > 0.3, and the galactic plane |b|< 10°. These cuts yielded a set of 66 610 stars, of which only 2440 had radial velocity measurements provided by Gaia.

To increase the amount of data with measured radial velocities and metallicities, we performed a series of data cross-matching steps. We used the SEGUE catalogue (Ahn et al. 2012) to complete the radial velocities and metallicities of our dataset, using the ones measured by Sloan (cross-match with ID table sdss_dr17.x1p5__specobj__gaia_dr3__gaia_source). In this case, no extra points were incorporated.

Additionally, we utilised the RR-Lyrae catalogue provided by Wang et al. (2022) to complete our catalogue with the sources not present in our dataset (288 points were incorporated) or to complete our radial velocity and metallicity data. After the cross-matching process, we ended up with a catalogue of 73 598 sources, of which 6523 have measured radial velocities. To obtain the radial velocity for the remaining sources, we applied a combined algorithm of data augmentation + random-forest regressor. In this case, the accuracy of our model is 49.0% (see Appendix B for details). Similar to the K-giant approach, we corrected the magnitudes to account for dust extinction. The absolute magnitude is connected to the metallicity through MG = 0.32[Fe/H] + 1.11 (Muraveva et al. 2018); therefore, one can obtain their distances using the distance modulus relationship.

To ensure consistency, we followed a similar approach to that used with the K-giants. To eliminate known objects from the Harris (1996) and Drlica-Wagner et al. (2020) catalogues, we applied an angular mask, taking into account the distance to the sources. Additionally, to exclude the Sagittarius stream from our analysis, we employed the cut-off criterion of |BSgr|< 15° for b > 0°. Then, we focused our analysis on objects with galactic distances between 30 kpc < Rgal < 100 kpc. In the final step, we removed stars that did not satisfy the condition μ α* 2 + ( μ δ +0.1) 2 < 0.5 2 $ \mu_{\alpha*}^2+(\mu_{\delta}+0.1)^2 < 0.5^2 $ (Conroy et al. 2021) after performing the solar reflex motion correction to the proper motions. Therefore, the final sample has 2446 sources.

Since our aim is to extract the mass of the DM subhalo surrounding the Large and Small Magellanic Clouds, we performed a transformation of coordinates to a new reference system. This particular coordinate reference system is centred on the centre of mass (CM) of the MCS, with the x axis aligned with the direction of the velocity of the CM (see Appendix A for details), and it is considered a rest frame. In order to obtain the position and velocity of the CM, we considered the LMC mass to be nine times the SMC mass (Craig et al. 2022).

3. Theoretical model and likelihood analysis

We used the theoretical model for stellar wakes from DM subhalos proposed by Buschmann et al. (2018), where the authors assumed a Plummer sphere for the density profile of the DM subhalo. The reference system used in this section corresponds to the one centred on the CM of the MCS, with the x axis in the direction of the velocity of the DM subhalo mass (see Appendix A.1 for details). From the collisionless Boltzmann equation, they derived the time-independent phase-space distribution function in the subhalo rest frame

f ( r ¯ , v ¯ , M s ) = f 0 ( v ¯ ) ( 1 + 2 G M s v 0 2 ( v ¯ + v ¯ s ) · α ¯ ) , α ¯ ( r ¯ , v ¯ , M s ) = 1 r v 1 + R s 2 r 2 1 + R s 2 r 2 v ¯ v r ¯ r 1 + R s 2 r 2 v ¯ · r ¯ rv . $$ \begin{aligned}&f\left(\bar{r},\bar{v},M_{\rm s}\right)=f_0(\bar{v}) \left(1+ \frac{2G M_{\rm s}}{v_0^2} \left(\bar{v}+\bar{v}_{\rm s}\right) \cdot \bar{\alpha }\right) , \\&\bar{\alpha }\left(\bar{r},\bar{v},M_{\rm s}\right)=\frac{1}{r v \sqrt{1+\frac{R_{\rm s}^2}{r^2}}} \frac{\sqrt{1+\frac{R_{\rm s}^2}{r^2}} \frac{\bar{v}}{v}-\frac{\bar{r}}{r}}{\sqrt{1+\frac{R_{\rm s}^2}{r^2}}-\frac{\bar{v}\cdot \bar{r}}{r v}} .\nonumber \end{aligned} $$(2)

In the equations, f 0 ( v ¯ ) = n 0 e ( v ¯ + v ¯ s ) 2 / v 0 2 π 3 / 2 v 0 3 $ f_0(\bar{v})=\frac{n_0 e^{-\left(\bar{v}+\bar{v}_{\mathrm{s}}\right)^2/v_0^2}}{\pi^{3/2} v_0^3} $, v 0 = 2 σ v $ {v_0= \sqrt{2} \sigma_{\mathrm{v}}} $ where σv is the velocity dispersion, v ¯ s $ \bar{v}_{\mathrm{s}} $ and Ms are the DM subhalo velocity and mass, respectively, R s ( M s ) = 1.62 M s / ( 10 8 M ) $ {R_{\mathrm{s}}(M_{\mathrm{s}})=1.62 \sqrt{M_{\mathrm{s}}/(10^8\,{M}_{\odot})}} $, and n0 is the star density inside the region of interest. This expression can be easily extended for three different velocity dispersions by including v0x (σvx), v0y (σvy), and v0z (σvz). In this case, the distribution function in the subhalo rest frame can be written as

f ( r ¯ , v ¯ , M s ) = n 0 e ( v x + v s v 0 x ) 2 ( v y v 0 y ) 2 ( v z v 0 z ) 2 π 3 / 2 v 0 x v 0 y v 0 z ( 1 + 2 G M s β ) , β ( r ¯ , v ¯ , M s ) = ( v x + v s v 0 x 2 i ̂ + v y v 0 y 2 j ̂ + v z v 0 z 2 k ̂ ) · α ¯ . $$ \begin{aligned}&f\left(\bar{r},\bar{v},M_{\rm s}\right)=\frac{n_0 \,e^{-\left(\frac{v_x+v_{\rm s}}{v_{0x}}\right)^2-\left(\frac{v_{ y}}{v_{0{ y}}}\right)^2-\left(\frac{v_z}{v_{0z}}\right)^2}}{\pi ^{3/2} v_{0x} v_{0{ y}} v_{0z}} \left(1+2G M_{\rm s} \beta \right) , \\&\beta \left(\bar{r},\bar{v},M_{\rm s}\right) = \left(\frac{v_x+v_{\rm s}}{v_{0x}^2} \hat{i} +\frac{v_{ y}}{v_{0{ y}}^2} \hat{j} +\frac{v_z}{v_{0z}^2} \hat{k} \right) \cdot \bar{\alpha }. \nonumber \end{aligned} $$(3)

In order to obtain the mass of the DM subhalo of the MCS, we performed a statistical analysis using the likelihood function to compare observational data in the new reference system and the theoretical model. This analysis was performed by using only the space data (3D) and the phase-space data (6D). The un-binned likelihood function for the 6D analysis is (Buschmann et al. 2018)

p 6 D ( M s , θ ) = e N s ( M s ) k = 1 N d f ( r ¯ k , v ¯ k , M s ) , $$ \begin{aligned} p_{6D}(M_{\rm s},\theta )= e^{-N_{\rm s}(M_{\rm s})} \prod _{k=1}^{N_{d}} f\left(\bar{r}_k,\bar{v}_k,M_{\rm s}\right), \end{aligned} $$(4)

where Nd is the number of stars in the region of interest (sphere of radius R centred on the CM of the MCS), Ns is the predicted number of stars in the same region, and θ are the fixed parameters of our model; that is, n0, v0, Rs, and v ¯ s $ \bar{v}_{\mathrm{s}} $. For the Plummer sphere and for the distribution function of Eq. (2), Ns(Ms) can be computed as (Buschmann et al. 2018)

N s ( M s ) = 4 3 π R 3 n 0 + 4 π G M s n 0 v 0 v s γ F ( v s v 0 ) , γ = R 2 1 + R s 2 R 2 R s 2 arcsinh ( R R s ) , F ( x ) = e x 2 0 x e y 2 d y . $$ \begin{aligned}&N_{\rm s}(M_{\rm s})=\frac{4}{3}\pi R^3 n_0 +\frac{4 \pi G M_{\rm s} n_0}{v_0 v_{\rm s}}\, \gamma \,F\left(\frac{v_{\rm s}}{v_0}\right) , \nonumber \\&\gamma =R^2\sqrt{1+\frac{R_{\rm s}^2}{R^2}} - R_{\rm s}^2\, \mathrm{arcsinh} \left(\frac{R}{R_{\rm s}}\right) ,\nonumber \\&F(x) = e^{-x^2} \int _0^x e^{{ y}^2} \mathrm{d}{ y} .\nonumber \end{aligned} $$

For the Plummer sphere and for the distribution function of Eq. (3), the predicted number of stars is

N s ( M s ) = 4 3 π R 3 n 0 + 4 G M s n 0 v 0 x v 0 y v 0 z γ I ( v ¯ 0 , v s ) , I ( v ¯ 0 , v s ) = d 3 p 2 π e p 2 cos ( 2 v s p x v 0 x ) ( ( p x v 0 x ) 2 + ( p y v 0 y ) 2 + ( p z v 0 z ) 2 ) 1 / 2 . $$ \begin{aligned}&N_{\rm s}(M_{\rm s}) = \frac{4}{3}\pi R^3 n_0 +\frac{4 G M_{\rm s} n_0}{v_{0x} v_{0{ y}} v_{0z}} \,\gamma \, I\left(\bar{v}_0,v_{\rm s}\right),\nonumber \\&I\left(\bar{v}_0,v_{\rm s}\right)=\int \frac{\mathrm{d}^3 p}{2\pi } \, e^{-p^2} \cos \left(\frac{2 v_{\rm s} p_x}{v_{0x}}\right) \nonumber \\&\quad \quad \quad \left(\left(\frac{p_x}{v_{0x}}\right)^2+\left(\frac{p_{ y}}{v_{0{ y}}}\right)^2+\left(\frac{p_z}{v_{0z}}\right)^2\right)^{-1/2}. \nonumber \end{aligned} $$

For the 3D data, the un-binned likelihood function is

p 3 D ( M s , θ ) = e N s ( M s ) k = 1 N d d 3 v f ( r ¯ k , v ¯ , M s ) . $$ \begin{aligned} p_{3D}(M_{\rm s},\theta )= e^{-N_{\rm s}(M_{\rm s})} \prod _{k=1}^{N_{\rm d}} \int \mathrm{d}^3 v \, f\left(\bar{r}_k,\bar{v},M_{\rm s}\right). \end{aligned} $$(5)

To determine the DM subhalo mass, we used a Markov-chain Monte Carlo (MCMC) method. For this purpose, we utilised the emcee package (Foreman-Mackey et al. 2013). We considered the function

λ ( M s ) = ln ( p x ( M s , θ ) ) , $$ \begin{aligned} \lambda (M_{\rm s})=\ln \left(p_{x}\left(M_{\rm s}, \theta \right)\right) , \end{aligned} $$(6)

where x stands for the 3D or 6D analysis. The prior used was 10 < log10Ms < 12. The functions, px, are presented in Eqs. (4) and (5). We considered 32 walkers and checked the convergence every 100 steps. To compute the uncertainties, we used the 16th, 50th, and 84th percentiles of the samples in the marginalised distributions (Foreman-Mackey et al. 2013).

4. Results

In Fig. 1, we present a Mollweide projection map displaying the distribution of our final dataset of 8504 stars in galactic coordinates (6058 K giants and 2446 RR Lyraes). To enhance the visual representation, the map has been smoothed using Gaussian functions with a full width at half maximum (FWHM) of 30°. The colour bar represents the density contrast, indicating the relative density variation from its mean value across the sample. The past 1 Gyr orbit of the CM of the MCS is shown with a magenta line, computed with the gala package using the MilkyWayPotential (Bovy 2015). The dotted blue and red lines represent the Pisces (Chandra et al. 2023b) and Virgo overdensities (Perottoni et al. 2022), respectively. The green lines are the polynomials adjusted by Chandra et al. (2023b) that limit the Pisces Plume. The black dots represent the members of the Magellanic Stellar Stream taken from Chandra et al. (2023a). The dark blue region corresponds to the masked region representing the galactic plane.

thumbnail Fig. 1.

Density distribution of K-giant and RR-Lyrae variables (Mollweide projection map) with 30 < Rgal < 60 kpc (top panel) and 60 < Rgal < 100 kpc (bottom panel). The data was smoothed using an FWHM of 30°. The magenta line represents the past orbit of the MCS CM. The dashed blue and red lines represent the Pisces and Virgo overdensities, respectively. The black dots represent the members of the Magellanic Stellar Stream. The green lines are the Pisces Plume.

Two distinct regions of overdensities can be observed. The first one, located in the northern hemisphere, with a longitude range between 225° and 315°, is associated with the collective response, also known as the global response. On the other hand, the southern feature appears to cover a larger area (−30° < l < 130°) and exhibits significant prominence at a longitude of 50° and a latitude of approximately −26°. This overdensity is associated with the local wake. A more minor component is also present in the northern hemisphere, within a longitude range of 30°–90°. It appears separate from the southern component due to the masking of the galactic plane, implemented to prevent contamination. The intensity of the wake is greater than that of the collective response. The ratio between the counts per pixel of the wake at l = 50°, b = −26° and the counts per pixel of the collective response at l = 279°, b = 24° (the coordinates of the highest overdensity of the collective response) is 1.29, considering the complete dataset. As one can see, the CM of the MCS past orbit is located over the local wake, and the deviations could be signalling the effect of the DM mass of the wake according to the results of Foote et al. (2023).

In the upper panel of Fig. 1, an inner region of the halo is shown (sources between 30 and 60 kpc from the galactic centre). It can be noted that both the wake and the collective do not correspond to either the Pisces or the Virgo overdensities. In previous works (Belokurov et al. 2019; Conroy et al. 2021), it was noticed that there exists a subregion produced by the Magellanic Clouds called the Pisces Plume. The southern overdensity that we identified as the wake does not fall within this region.

On the other hand, in the lower panel, we plot the outer region of the halo, between 60 and 100 kpc, a region studied by Belokurov et al. (2019) and Conroy et al. (2021). Once again, the collective does not belong to the Virgo overdensity. Nevertheless, the global response could be truncated due to the masking of the galactic plane and the Sagittarius stream. However, part of the wake lies on the edge of the Pisces overdensity, and the maximum of the wake is indeed located in the Pisces Plume. Additionally, it has been verified that none of the catalogued points belonging to the Magellanic Stellar Stream (Chandra et al. 2023a) are found in our dataset.

Comparing our results with the ones obtained by Conroy et al. (2021), we observe some slight differences in the locations of the overdensities. Specifically, our sample’s maximum southern overdensity is slightly displaced further north. Similarly, the maximum northern overdensity in our sample is slightly shifted towards the south and east. In particular, we have also compared only our K-giant sample results with the final public catalogue developed by Conroy et al. (2021). The coordinates of the maximum overdensities of the local wake and the collective response for our dataset are (l = 52°, b = −24°) and (l = 275°, b = 25°), respectively; meanwhile, for Conroy’s data they are approximately (l = 49°, b = −54°) for the wake and (l = 326°, b = 54°) for the collective. Our definition of the local wake is larger than the one proposed by Conroy et al. (2021), and the ratio between the counts per pixel of the wake and the counts per pixel of the collective response at each maximum is 1.5 for our data and 1.33 for Conroy’s data-set. However, if we consider only K giants located within 60 < Rgal < 100 kpc, we successfully replicated the wake’s position, with its peak occurring at l = 57°, b = −51°. Moreover, when we compare our map with the simulations presented in Conroy et al. (2021), we observe quite an agreement regarding the positions of the overdensities.

In Fig. 2, we plot the superficial overdensity of the wake along with the past orbit of the CM (magenta line), in a new coordinate frame; namely, the orbit frame (see Appendix A.2). In this frame, the plane x* − y* contains the CM orbit and the z* axis is perpendicular to the CM orbit. The origin of this new coordinate frame is the current location of the CM, and the x* axis is coincident with the direction of the DM subhalo mass velocity. As one can see, the wake is located in x* < 0 and it moves towards the perturber. In particular, the lower panel of Fig. 2 is in good agreement with the results presented in Fig. 1 of Buschmann et al. (2018).

thumbnail Fig. 2.

Overdensity as a function of the position. Magenta line: past orbit of the CM. Magenta star: current position of the CM. Pink arrow: velocity of the CM. Red arrow: mean velocity of the wake or collective. Green arrow: mean velocity of a bin of 15 kpc. Cyan star: position of the galactic centre. Top panel: orbit plane with |z*|< 10 kpc. Bottom panel: x* − z* plane (perpendicular to the orbit plane).

The maximum value of the density is approximately x* = −48.97 kpc, y* = 55.65 kpc, and z* = −7.94 kpc in the new orbit frame. Using the stars in a 10 kpc neighbourhood (177 stars), we computed the velocity dispersion, resulting in (48 ± 3) km s−1. The characterisation of its complete stellar population and dynamical properties will be addressed in a forthcoming work.

The Gaia satellite data and spectroscopic surveys have unveiled the structure, composition, and formation history of the MW’s stellar halo in recent years (Helmi et al. 2018; Belokurov et al. 2019; Kruijssen et al. 2020; Callingham et al. 2022). Satellite galaxies, globular clusters, stars, and streams are now associated with different halo components, which constitute the remains of our galaxy’s past merger events.

In order to study the possible origin of the stellar populations of the wake and collective response, we have made an E vs. Lz diagram (see Fig. 3). We have taken the same coordinate system convention as the one adopted by Callingham et al. (2022). In this diagram, we have discriminated the different known mergers of the MW (colour points, data extracted from Callingham et al. 2022), the LMC (star), the SMC (squared), and the wake and collective response (cyan and magenta regions, respectively). As one can see, both the wake and the collective are extended regions in the diagram without a defined sign of Lz but limited in energy; however, it is not in the range of the Gaia-Sausage-Enceladus energy. The mean values for the wake are E = −0.703 × 105 km2s−2 and Lz = 890.433 km kpc s−1; meanwhile, for the collective, they are E = −0.716 × 105 km2 s−2 and Lz = 902.391 km kpc s−1.

thumbnail Fig. 3.

E − Lz diagram. The mergers data were taken from Callingham et al. (2022). Cyan region: wake. Magenta region: collective response.

It should be noted that the stellar population of the MW’s Halo was accreted by mergers of older origin such as Sequoia, Saggitarius, Helmi, ED-3-4-5-6, Typhon, and L-RL64 (Dodd et al. 2023). Consequently, the recent impact of the LMC could dynamically affect all these stellar populations; however, more studies are needed.

Next, we performed a statistical analysis in order to obtain the DM subhalo mass of the CM, and therefore the DM subhalo mass of the LMC. The radius for the region of interest was fixed at R = 100 kpc, the subhalo velocity, vs, was fixed at 314.23 km s−1 (van der Marel et al. 2002; Martínez-Delgado et al. 2019), and n0 was obtained from the reduced observational data described in the previous section. For the velocity dispersion, we performed a statistical analysis of the data and obtained the velocity standard deviation in each axis. For the analysis using Eq. (2), we considered the larger component of v ¯ 0 $ \bar{v}_0 $ in the calculations. We used the described density profile in Sect. 3 and we present our results in Fig. 4, along with the DM subhalo mass estimation of the LMC published in the literature (Watkins et al. 2024; Koposov et al. 2023; Shipp et al. 2021; Vasiliev et al. 2021; Erkal et al. 2019; Peñarrubia et al. 2016; using a different method, indicated with colours) along with our results (the last three values, with their corresponding statistical errors).

thumbnail Fig. 4.

Dark matter subhalo mass estimation of the LMC. Orange line: kinematic estimation from MW satellites. Red lines: estimation from stellar streams. Green line: estimation based on momentum balance in the Local Group. Blue lines: our results obtained from the likelihood analysis. Case (a): 3D data, M LMC = 2 . 289 0.240 + 0.260 × 10 11 M $ M_{\mathrm{LMC}}=2.289^{+0.260}_{-0.240}\times 10^{11}\,{M}_{\odot} $. Case (b): 6D data (Eq. (3)), M LMC = 1 . 787 0.069 + 0.072 × 10 11 M $ M_{\mathrm{LMC}}=1.787^{+0.072}_{-0.069}\times 10^{11}\,{M}_{\odot} $. Case (c): 6D data (Eq. (2)), M LMC = 1 . 686 0.072 + 0.071 × 10 11 M $ M_{\mathrm{LMC}}=1.686^{+0.071}_{-0.072}\times 10^{11}\,{M}_{\odot} $. The error bars of our results are purely statistical and based on a restricted one-parameter model.

The DM subhalo mass was computed by using the space distribution function (case (a)) and the phase-space distribution function of Eq. (3) (case (b)) and Eq. (2) (case (c)). Our results are consistent, despite the distribution function used in the analysis. However, the fit obtained using only the space data is slightly higher than the 6D analysis results. Furthermore, our findings agree within 3σ with the literature (Vasiliev 2023).

We also performed the statistical analysis (6D Eq. (2)) using only the data with measured radial velocity (Gaia Collaboration 2023; Ahn et al. 2012; Wang et al. 2022). We found an LMC subhalo mass of M LMC = 1 . 594 0.196 + 0.203 × 10 11 M $ M_{\mathrm{LMC}}=1.594^{+0.203}_{-0.196}\times 10^{11}\,{M}_{\odot} $, which is in good agreement with our results.

5. Conclusions

In this work, we employed the recently published Gaia Data Release 3, which improves the precision of proper motions along with the Segue catalogue (Ahn et al. 2012) and the one provided by Wang et al. (2022). This enabled us to extend the K-giant catalogue originally provided by Conroy et al. (2021), and also to construct a catalogue for RR-Lyrae stars, both in 6D data. We reproduced the previously published results and identified the overdensities associated with a wake and the collective response using these two halo tracers. A notable finding of this study is the extension of the southern overdensity towards lower galactocentric distances; that is, between 30 and 100 kpc. Moreover, we were able to show that the southern overdensity, identified as the wake, trails the CM of Magellanic Clouds (see Fig. 2).

We have confirmed that the Pisces plume overdensity, described in Belokurov et al. (2019), is associated with the wake of the Magellanic Clouds in the outer regions of the MW’s halo. Furthermore, we have discovered that the overdensity on the halo’s stellar population, caused by the Magellanic Clouds’ wake and the global response, affects the stars in the MW’s halo, regardless of which past merger event they were accreted from.

As per the theoretical proposal made by Buschmann et al. (2018), we were able to estimate the mass of the LMC DM subhalo for the first time by using Gaia observational data. We found a reliable estimation of the DM halo surrounding the LMC by performing two different analyses, using only the space distribution data and using both the phase and space data. Considering a relationship between the Large and Small Magellanic Clouds’ masses, our study has successfully determined the mass of the DM subhalo of the larger cloud. Even more, our findings are in agreement with prior results (Correa Magnus & Vasiliev 2021; Koposov et al. 2023; Shipp et al. 2021; Vasiliev et al. 2021; Erkal et al. 2019; Peñarrubia et al. 2016), within 3σ. This consistency with previous studies indicates the reliability of our methodology. Additionally, this method gives competitive errors compared to different mass determination methods. It is important to point out that the errors mentioned were calculated using the assumptions of a Plummer profile for a spherical subhalo. However, if a more complex multi-parameter model such as an ellipsoidal Navarrro, Frenk and White (NFW) model were used, it is anticipated that the errors would increase.

Acknowledgments

This work was partially supported by grants from the National Research Council of Argentina (CONICET PIP 616). K.J.F. is a Post Doctoral fellow of the CONICET. M.E.M. and M.D. are members of the Scientific Research Career of the CONICET. M.D. work was supported by the Preparing for Astrophysics with LSST Program, funded by the Heising Simons Foundation through grant 2021-2975, and administered by Las Cumbres Observatory. We would like to express our gratitude to the referee for providing us with their valuable comments and suggestions. M.D. thanks valuable suggestions by Nicolás Garavito-Camargo and the support of the CCA and the Flatiron Institute. M.D. also thanks the possibility of virtual participation in the following programs: “Building a physical understanding of galaxy evolution with data-driven astronomy” at the KITP-UCB, the Chicago “Gaia DR3 Sprint” and “Streams22: Community Atlas of Tidal Streams.” at the Carnegie Observatories and the Flatiron Institute. K.J.F wants to thank Alejo Molina Lera and Martín Gamboa Lerena for the useful discussions. 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), and code developed by the Gaia Project Scientist Support Team. Funding for the DPAC has been provided by national institutions, in particular, the institutions participating in the Gaia Multilateral Agreement. This research has utilized the following software: Astropy (Astropy Collaboration 2013), Matplotlib (Hunter 2007), Pandas (McKinney et al. 2010), Seaborn (Waskom 2021), Healpy (Zonca et al. 2019), SciPy (Virtanen et al. 2020), NumPy (Harris et al. 2020), emcee (Foreman-Mackey et al. 2013), gala (Price-Whelan 2017; Price-Whelan et al. 2020), The Jupyter Notebook (Kluyver et al. 2016), Scikit-learn (Pedregosa et al. 2011), Pzflow (Crenshaw et al. 2024), and TOPCAT (Taylor 2005).

References

  1. Abdalla, H., Aharonian, F., Benkhali, F. A., et al. 2022, Phys. Rev. Lett., 129, 111101 [NASA ADS] [CrossRef] [Google Scholar]
  2. Abe, H., Abe, S., Acciari, V. A., et al. 2023, Phys. Rev. Lett., 130, 061002 [NASA ADS] [CrossRef] [Google Scholar]
  3. Acharyya, A., Archer, A., Bangale, P., et al. 2023, ApJ, 945, 101 [NASA ADS] [CrossRef] [Google Scholar]
  4. Aguilar-Argüello, G., Valenzuela, O., & Trelles, A. 2022, A&A, 663, A93 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Ahn, C. P., Alexandroff, R., Allende Prieto, C., et al. 2012, ApJS, 203, 21 [Google Scholar]
  6. Amaré, J., Cebrián, S., Cintas, D., et al. 2022, Moscow Univ. Phys. Bull., 77, 322 [CrossRef] [Google Scholar]
  7. Aprile, E., Abe, K., Agostini, F., et al. 2023, Phys. Rev. Lett., 131, 041003 [CrossRef] [PubMed] [Google Scholar]
  8. Astropy Collaboration (Robitaille, T. P., et al.) 2013, A&A, 558, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Barberio, E., Baroncelli, T., Bignell, L. J., et al. 2023, Eur. Phys. J. C, 83, 878 [NASA ADS] [CrossRef] [Google Scholar]
  10. Bazarov, A., Benito, M., Hütsi, G., et al. 2022, Astron. Comput., 41, 100667 [NASA ADS] [CrossRef] [Google Scholar]
  11. Belokurov, V., Deason, A. J., Erkal, D., et al. 2019, MNRAS, 488, L47 [Google Scholar]
  12. Bennett, M., & Bovy, J. 2019, MNRAS, 482, 1417 [NASA ADS] [CrossRef] [Google Scholar]
  13. Bernabei, R., Belli, P., Bussolotti, A., et al. 2022, in The Fifteenth Marcel Grossmann Meeting on General Relativity, eds. E. S. Battistelli, R. T. Jantzen, & R. Ruffini, 1285 [CrossRef] [Google Scholar]
  14. Bovy, J. 2015, ApJS, 216, 29 [NASA ADS] [CrossRef] [Google Scholar]
  15. Breiman, L. 2001, Mach. Learn., 45, 5 [Google Scholar]
  16. Buehler, R., & Desjacques, V. 2023, Phys. Rev. D, 107, 023516P [NASA ADS] [CrossRef] [Google Scholar]
  17. Buschmann, M., Kopp, J., Safdi, B. R., & Wu, C.-L. 2018, Phys. Rev. Lett., 120, 211101 [NASA ADS] [CrossRef] [Google Scholar]
  18. Callingham, T. M., Cautun, M., Deason, A. J., et al. 2022, MNRAS, 513, 4107 [NASA ADS] [CrossRef] [Google Scholar]
  19. Carr, B., & Kühnel, F. 2020, Ann. Rev. Nucl. Part. Scie., 70, 355 [NASA ADS] [CrossRef] [Google Scholar]
  20. Chandra, V., Naidu, R. P., Conroy, C., et al. 2023a, ApJ, 956, 110 [CrossRef] [Google Scholar]
  21. Chandra, V., Naidu, R. P., Conroy, C., et al. 2023b, ApJ, 951, 26 [NASA ADS] [CrossRef] [Google Scholar]
  22. Choi, J., Dotter, A., Conroy, C., et al. 2016, ApJ, 823, 102 [Google Scholar]
  23. Clowe, D., Bradac, M., Gonzalez, A. H., et al. 2006, ApJ, 648, L109 [NASA ADS] [CrossRef] [Google Scholar]
  24. Conroy, C., Naidu, R. P., Garavito-Camargo, N., et al. 2021, Nature, 592, 534 [NASA ADS] [CrossRef] [Google Scholar]
  25. Correa Magnus, L., & Vasiliev, E. 2021, MNRAS, 511, 2610 [Google Scholar]
  26. Craig, P. A., Chakrabarti, S., Baum, S., & Lewis, B. T. 2022, MNRAS, 517, 1737 [NASA ADS] [CrossRef] [Google Scholar]
  27. Crenshaw, J. F., Connolly, A., & Kalmbach, B. 2021, Am. Astron. Soc. Meet. Abstr., 53, 230.01 [Google Scholar]
  28. Crenshaw, J. F., Connolly, A., & Kalmbach, B. 2024, https://doi.org/10.5281/zenodo.10710271 [Google Scholar]
  29. Cunningham, E. C., Garavito-Camargo, N., Deason, A. J., et al. 2020, ApJ, 898, 4 [CrossRef] [Google Scholar]
  30. Dodd, E., Callingham, T. M., Helmi, A., et al. 2023, A&A, 670, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Dotter, A. 2016, ApJS, 222, 8 [Google Scholar]
  32. Drimmel, R., & Poggio, E. 2018, Res. Notes Am. Astron. Soc., 2, 210 [Google Scholar]
  33. Drlica-Wagner, A., Bechtol, K., Mau, S., et al. 2020, ApJ, 893, 47 [CrossRef] [Google Scholar]
  34. Durkan, C., Bekasov, A., Murray, I., & Papamakarios, G. 2019, arXiv e-prints [arXiv:1906.04032] [Google Scholar]
  35. Erkal, D., Belokurov, V., Laporte, C. F. P., et al. 2019, MNRAS, 487, 2685 [Google Scholar]
  36. Foote, H. R., Besla, G., Mocz, P., et al. 2023, ApJ, 954, 163 [CrossRef] [Google Scholar]
  37. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
  38. Furlanetto, S. R., & Loeb, A. 2002, ApJ, 565, 854 [CrossRef] [Google Scholar]
  39. Gaia Collaboration (Prusti, T., et al.) 2016, A&A, 595, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  40. Gaia Collaboration (Vallenari, A., et al.) 2023, A&A, 674, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Garavito-Camargo, N., Besla, G., Laporte, C. F. P., et al. 2019, ApJ, 884, 51 [NASA ADS] [CrossRef] [Google Scholar]
  42. GRAVITY Collaboration (Abuter, R., et al.) 2019, A&A, 625, L10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Green, G. 2018, J. Open Source Softw., 3, 695 [NASA ADS] [CrossRef] [Google Scholar]
  44. Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [NASA ADS] [CrossRef] [Google Scholar]
  45. Harris, W. E. 1996, AJ, 112, 1487 [Google Scholar]
  46. Helmi, A., Babusiaux, C., Koppelman, H. H., et al. 2018, Nature, 563, 85 [Google Scholar]
  47. Hui, L., Ostriker, J. P., Tremaine, S., & Witten, E. 2017, Phys. Rev. D, 95, 043541 [NASA ADS] [CrossRef] [Google Scholar]
  48. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [NASA ADS] [CrossRef] [Google Scholar]
  49. Katz, D., Sartoretti, P., Guerrier, A., et al. 2023, A&A, 674, A5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  50. Kluyver, T., Ragan-Kelley, B., Pérez, F., et al. 2016, Jupyter Notebooks – a Publishing Format for Reproducible Computational Workflows (IOS Press), 87 [Google Scholar]
  51. Koposov, S. E., Erkal, D., Li, T. S., et al. 2023, MNRAS, 521, 4936 [NASA ADS] [CrossRef] [Google Scholar]
  52. Kruijssen, J. M. D., Pfeffer, J. L., Chevance, M., et al. 2020, MNRAS, 498, 2472 [NASA ADS] [CrossRef] [Google Scholar]
  53. Martínez-Delgado, D., Vivas, A. K., Grebel, E. K., et al. 2019, A&A, 631, A98 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Massey, R., Kitching, T., & Richard, J. 2010, Rep. Progr. Phys., 73, 086901 [CrossRef] [Google Scholar]
  55. McConnachie, A. W. 2012, AJ, 144, 4 [Google Scholar]
  56. McKinney, W. 2010, in Proceedings of the 9th Python in Science Conference, eds. S. van der Walt, & J. Millman, 51 [Google Scholar]
  57. Mo, H., van den Bosch, F., & White, S. 2010, Galaxy Formation and Evolution (Cambridge: Cambridge University Press) [CrossRef] [Google Scholar]
  58. Muraveva, T., Delgado, H. E., Clementini, G., Sarro, L. M., & Garofalo, A. 2018, MNRAS, 481, 1195 [Google Scholar]
  59. Naik, A. P., & Widmark, A. 2023, MNRAS, 527, 11559 [NASA ADS] [CrossRef] [Google Scholar]
  60. Paxton, B., Bildsten, L., Dotter, A., et al. 2011, ApJS, 192, 3 [Google Scholar]
  61. Pedregosa, F., Varoquaux, G., Gramfort, A., et al. 2011, J. Mach. Learn. Res., 12, 2825 [Google Scholar]
  62. Peñarrubia, J., Gómez, F. A., Besla, G., Erkal, D., & Ma, Y.-Z. 2016, MNRAS, 456, L54 [CrossRef] [Google Scholar]
  63. Perottoni, H. D., Limberg, G., Amarante, J. A. S., et al. 2022, ApJ, 936, L2 [NASA ADS] [CrossRef] [Google Scholar]
  64. Planck Collaboration VI. 2020, A&A, 641, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Price-Whelan, A. M. 2017, J. Open Source Softw., 2, 388 [NASA ADS] [CrossRef] [Google Scholar]
  66. Price-Whelan, A., Sipőcz, B., Lenz, D., et al. 2020, https://doi.org/10.5281/zenodo.4159870 [Google Scholar]
  67. Riello, M., De Angeli, F., Evans, D. W., et al. 2021, A&A, 649, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  68. Rubin, V. C., Ford, W., & Kent, J. 1970, ApJ, 159, 379 [NASA ADS] [CrossRef] [Google Scholar]
  69. Shipp, N., Erkal, D., Drlica-Wagner, A., et al. 2021, ApJ, 923, 149 [NASA ADS] [CrossRef] [Google Scholar]
  70. Taylor, M. B. 2005, in Astronomical Data Analysis Software and Systems XIV, eds. P. Shopbell, M. Britton, & R. Ebert, ASP Conf. Ser., 347, 29 [Google Scholar]
  71. Tulin, S., & Yu, H.-B. 2018, Phys. Rep., 730, 1 [Google Scholar]
  72. van der Marel, R. P., Alves, D. R., Hardy, E., & Suntzeff, N. B. 2002, AJ, 124, 2639 [NASA ADS] [CrossRef] [Google Scholar]
  73. Vasiliev, E. 2023, Galaxies, 11, 59 [NASA ADS] [CrossRef] [Google Scholar]
  74. Vasiliev, E., Belokurov, V., & Erkal, D. 2021, MNRAS, 501, 2279 [NASA ADS] [CrossRef] [Google Scholar]
  75. Villanueva-Domingo, P., Mena, O., & Palomares-Ruiz, S. 2021, Front. Astron. Space Sci., 8 [Google Scholar]
  76. Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Methods, 17, 261 [Google Scholar]
  77. Wang, F., Zhang, H.-W., Xue, X.-X., et al. 2022, MNRAS, 513, 1958 [NASA ADS] [CrossRef] [Google Scholar]
  78. Waskom, M. L. 2021, J. Open Source Softw., 6, 3021 [CrossRef] [Google Scholar]
  79. Watkins, L. L., van der Marel, R. P., & Bennet, P. 2024, ApJ, 963, 84 [NASA ADS] [CrossRef] [Google Scholar]
  80. Weinberg, M. D. 1986, ApJ, 300, 93 [NASA ADS] [CrossRef] [Google Scholar]
  81. Zavala, J., & Frenk, C. S. 2019, Galaxies, 7, 81 [NASA ADS] [CrossRef] [Google Scholar]
  82. Zonca, A., Singer, L., Lenz, D., et al. 2019, J. Open Source Softw., 4, 1298 [Google Scholar]
  83. Zwicky, F. 1933, Helv. Phys. Acta, 6, 110 [NASA ADS] [Google Scholar]

Appendix A: Large Magellanic Cloud rest frame

A.1. Centre-of-mass rest frame

The coordinate system used to compare the data with the theoretical model has its origin in the LMC and SMC CM and the x axis orientated according to the DM subhalo velocity, v ¯ s $ \bar{v}_{\mathrm{s}} $. In order to obtain the coordinates of our dataset in such a rest frame, we performed the following steps:

  1. we boosted the data to the new frame by r ¯ boost = r ¯ obs r ¯ cm $ {\bar{r}_{boost}=\bar{r}_{obs}-\bar{r}_{cm}} $;

  2. we performed the rotations upon the boosted data using the matrix

M = ( cos θ 1 cos θ 2 sin θ 1 cos θ 2 sin θ 2 sin θ 1 cos θ 1 0 cos θ 1 sin θ 2 sin θ 1 sin θ 2 cos θ 2 ) , $$ \begin{aligned} M&=\left( \begin{array}{ccc} \cos \theta _1 \cos \theta _2&\sin \theta _1 \cos \theta _2&\sin \theta _2 \\ -\sin \theta _1&\cos \theta _1&0 \\ - \cos \theta _1 \sin \theta _2&-\sin \theta _1 \sin \theta _2&\cos \theta _2 \end{array} \right)\,, \end{aligned} $$(A.1)

to obtain the coordinate (x′, y′, z′) in the new rest frame. The velocity had to be transformed as well and the final velocity had to be boosted using v ¯ s = v s i ̂ $ \bar{v}_{\mathrm{s}}=v_{\mathrm{s}} \hat{i}^{\prime} $. The angles are defined as

tan θ 1 = ( v cm ) y ( v cm ) x , $$ \begin{aligned} \tan \theta _1&=\frac{\left(v_{cm}\right)_y}{\left(v_{cm}\right)_x} \,, \end{aligned} $$(A.2)

tan θ 2 = ( v cm ) z ( v cm ) x 2 + ( v cm ) y 2 , $$ \begin{aligned} \tan \theta _2&=\frac{\left(v_{cm}\right)_z}{\sqrt{\left(v_{cm}\right)_x^2+\left(v_{cm}\right)_y^2}}\,, \end{aligned} $$(A.3)

where ((vcm)x, (vcm)y, (vcm)z) is the CM velocity in the solar coordinate system.

A.2. Orbit frame

The coordinate system used to compare our results with the Fig. 1 presented by Buschmann et al. (2018) has its origin in the CM of the MCS, its z* axis perpendicular to the orbital plane, and the x* axis orientated in the direction of the velocity of the DM subhalo. To obtain these new coordinates from the CM rest frame, we had to perform a rotation of an angle, θorb, according to the following matrix:

M orb = ( 1 0 0 0 cos θ orb sin θ orb 0 sin θ orb cos θ orb ) , $$ \begin{aligned} M_{orb}&= \left( \begin{array}{ccc} 1&0&0 \\ 0&\cos \theta _{orb}&\sin \theta _{orb} \\ 0&-\sin \theta _{orb}&\cos \theta _{orb} \end{array} \right)\,, \end{aligned} $$(A.4)

where tan θ orb = z ¯ orb / y ¯ orb $ \tan \theta_{orb}=\overline{z}_{orb}/\overline{y}_{orb} $, z ¯ orb $ \overline{z}_{orb} $ ( y ¯ orb ) $ \left(\overline{y}_{orb}\right) $ is the mean value of the z (y) coordinate in the CM rest frame of the CM’s past orbit.

Appendix B: Estimation of radial velocities using machine learning

To complete the phase information for all the halo stars in our sample, we applied two machine learning techniques to assign radial velocities to those stars without such measurements.

B.1. K-giant radial velocity

Our first sample consists of 245086 K-giant stars, of which 10989 have measured radial velocities. The last subsample was used to train a random forest regressor (RF) (Breiman 2001). The chosen predictors for this study were the angular coordinates, proper motions, G magnitudes, BP and RP colours, distances to the Sun and the galactic centre, and the Galactocentric Cartesian coordinates. To prevent overfitting, a standard cross-validation analysis was performed. The RF hyper-parameters were tuned using GridSearchCV from the Scikit Learn library, resulting in the following values: [n_estimators = 4900, random_state = 0]

B.2. RR-Lyrae radial velocity

We observed a decrease in the number of measured radial velocities for RR-Lyrae stars, resulting in a drop in the fraction of measured radial velocities to (5510/67276), considering galactocentric distances between 10 to 100 kpc. To tackle this issue, we aimed to model the spatial distribution of radial velocities. To achieve this, we utilised normalising flows (NFs) (Durkan et al. 2019), as have been implemented by Crenshaw et al. (2021), to model the joint posterior probability of radial velocities and predictors using a subsample of the features described earlier. However, we limited the feature space to prevent any bias against less luminous stars. Therefore, we did not consider the magnitude and colours of stars as predictors in this case.

We used probabilistic modelling to generate a radial velocity distribution outcome for a significant number of stars (100000). Therefore, the normalizing flow was used to augment and generalise the training dataset. We evaluated the marginal probability of the radial velocity given the values of other variables (predictors) and obtained a vector of probabilities. To complement this method, we used the RF algorithm to map the posterior conditional radial velocity distribution to the real measured value in the training sample, similar to the method used for K-giant stars. The RF hyper-parameters were tuned using GridSearchCV from the Scikit Learn library, resulting in the following values: [max_depth = 50, max_features = 8, min_samples_leaf = 1, min_samples_split = 6, n_estimators = 800]. It is important to note that the result of the NFs is a vector that represents the conditional radial velocity. The vector is measured on a 1000-dimensional grid that samples velocities ranging from -700 to 700 km/s. To prevent overfitting, cross-validation was employed by splitting the data in two and using 80% for training and 20% for validation, similar to the previous case.

B.3. Goodness of fit in radial velocity regression

R-Squared (R2), or the coefficient of determination, is a statistical measure used to determine the proportion of variance in the dependent variable that can be explained by the independent variable in a regression model. The statistics were calculated for our two samples of stars, resulting in values of 0.86 for K-giant stars and 0.5 for RR-Lyrae ones, respectively. These values are comparable to those recently reported by Naik & Widmark (2023) using Bayesian neural networks, and provide us with a complete sample of halo stars in the phase space.

We tested the inferred stars’ radial velocities by measuring the mean radial velocities of satellite galaxies and compared them with the tabulated velocity values (McConnachie 2012). The results are presented in Table. B.1. Our machine-learning results show that the average velocity of the stars on these satellite galaxies is well reproduced. The relevant error of the average velocity for the Sculptor galaxy is due to the lack of any measured radial velocity on this object.

Table B.1.

Comparison between the mean radial velocity inferred with machine learning and the values tabulated for MW satellite galaxies (McConnachie 2012).

Using this method, we can estimate the mean velocity field in the halo, which is necessary for the likelihood estimation of the subhalo mass of the LMC.

Appendix C: Validation of estimated distances for K-giant stars

To test the photometric distance obtained for the K giants, we computed the photometric distances of the globular clusters NGC7006, NGC5694, NGC2419, and NGC6229, and of the LMC, SMC, Carina, Draco, and Sculptor.

The observational data was extracted from Gaia Data Release 3 (Gaia Collaboration 2023, 2016), using the option single object searcher to avoid field-contamination. Additionally, the data for the LMC and the SMC were selected on a circular disk in (l, b) centred on the location of each MC, with a radius of 2°. For Carina, Draco, and Sculptor we selected data using an angular mask of three times the tidal radius reported in Drlica-Wagner et al. (2020). We performed the data reduction or treatment indicated in the text for the K giant (that is, corrected for dust extinction, discarded all sources with E(B − V) > 0.3, corrected the magnitudes, and performed the 3σ cut upon the corrected BP and RP flux excess factor (C*)). Finally, we selected the giant branch and performed the cross-match with the K-giant catalogue given by Gaia (gaiadr3.astrophysical parameters). For the LMC and SMC, we also restricted the dataset in the reported Gaia parallax.

We calculated the photometric distance of each source and the mean value of each object. The results are shown in Fig. C.1 and, as can be seen, the computed photometric distances and the values reported in the literature (Harris 1996; Drlica-Wagner et al. 2020) are in good agreement.

thumbnail Fig. C.1.

Comparison between the photometric distances computed using Eq. (1) and the tabulated distances (Harris 1996; Drlica-Wagner et al. 2020). Solid olive line: photometric distance equals tabulated distance. Dashed brown line: ±10% with respect to the one-to-one line. Solid purple line: ±20% with respect to the one-to-one line. The vertical lines are the statistical errors.

We have also applied our fit to K-giant stars from the Conroy et al. (2021) catalogue and found that our calculated distances are in agreement with the distances reported by Conroy et al. (2021) within 10%.

All Tables

Table B.1.

Comparison between the mean radial velocity inferred with machine learning and the values tabulated for MW satellite galaxies (McConnachie 2012).

All Figures

thumbnail Fig. 1.

Density distribution of K-giant and RR-Lyrae variables (Mollweide projection map) with 30 < Rgal < 60 kpc (top panel) and 60 < Rgal < 100 kpc (bottom panel). The data was smoothed using an FWHM of 30°. The magenta line represents the past orbit of the MCS CM. The dashed blue and red lines represent the Pisces and Virgo overdensities, respectively. The black dots represent the members of the Magellanic Stellar Stream. The green lines are the Pisces Plume.

In the text
thumbnail Fig. 2.

Overdensity as a function of the position. Magenta line: past orbit of the CM. Magenta star: current position of the CM. Pink arrow: velocity of the CM. Red arrow: mean velocity of the wake or collective. Green arrow: mean velocity of a bin of 15 kpc. Cyan star: position of the galactic centre. Top panel: orbit plane with |z*|< 10 kpc. Bottom panel: x* − z* plane (perpendicular to the orbit plane).

In the text
thumbnail Fig. 3.

E − Lz diagram. The mergers data were taken from Callingham et al. (2022). Cyan region: wake. Magenta region: collective response.

In the text
thumbnail Fig. 4.

Dark matter subhalo mass estimation of the LMC. Orange line: kinematic estimation from MW satellites. Red lines: estimation from stellar streams. Green line: estimation based on momentum balance in the Local Group. Blue lines: our results obtained from the likelihood analysis. Case (a): 3D data, M LMC = 2 . 289 0.240 + 0.260 × 10 11 M $ M_{\mathrm{LMC}}=2.289^{+0.260}_{-0.240}\times 10^{11}\,{M}_{\odot} $. Case (b): 6D data (Eq. (3)), M LMC = 1 . 787 0.069 + 0.072 × 10 11 M $ M_{\mathrm{LMC}}=1.787^{+0.072}_{-0.069}\times 10^{11}\,{M}_{\odot} $. Case (c): 6D data (Eq. (2)), M LMC = 1 . 686 0.072 + 0.071 × 10 11 M $ M_{\mathrm{LMC}}=1.686^{+0.071}_{-0.072}\times 10^{11}\,{M}_{\odot} $. The error bars of our results are purely statistical and based on a restricted one-parameter model.

In the text
thumbnail Fig. C.1.

Comparison between the photometric distances computed using Eq. (1) and the tabulated distances (Harris 1996; Drlica-Wagner et al. 2020). Solid olive line: photometric distance equals tabulated distance. Dashed brown line: ±10% with respect to the one-to-one line. Solid purple line: ±20% with respect to the one-to-one line. The vertical lines are the statistical errors.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.