Discovery of NES, an extended tidal structure in the North-East of the Large Magellanic Cloud

We report on the discovery of a new diffuse stellar sub-structure protruding for 5 degrees from the North-Eastern rim of the LMC disc. The structure, that we dub North-East Structure (NES), was identified by applying a Gaussian Mixture Model to a sample of strictly selected candidate members of the Magellanic System, extracted from the Gaia EDR3 catalogue. The NES fills the gap between the outer LMC disk and other known structures in the same region of the LMC, namely the Northern tidal arm (NTA) and the Eastern sub-structures (ES). Particularly noteworthy is that the NES is placed in a region where N-body simulations foresee a bending of the LMC disc due to tidal stresses induced by the MW. The velocity field in the plane of the sky indicates that the complex of tidal structures in the North-Eastern part of the LMC, including NES, is subject to coherent radial motions. Additional data, as well as extensive dynamical modeling, is required to shed light on the origin of NES as well as on the relationships with the surrounding substructures.


INTRODUCTION
The Large Magellanic Cloud (LMC) is the largest dwarf satellite of the Milky Way (MW) and along with its smaller companion, the Small Magellanic Cloud (SMC), represents a very close example of a physical system in an ongoing tidal stripping process. Some signatures of their mutual interaction are very prominent in the sky, such as the Magellanic stream, a cloud of HI gas that extends for more than 180 • in the Southern hemisphere (Putman et al. 2003;Brüns et al. 2005;D'Onghia & Fox 2016), or its counterpart, the Leading arm (Putman et al. 1998). The Magellanic Clouds (MCs) are also connected by the Magellanic bridge, a stream likely originated as a consequence of the last close encounter between the MCs (Zivick et al. 2018), made of HI gas (Hindman et al. 1963), young (e.g., Harris 2007;Skowron et al. 2014;Noël et al. 2015;Mackey et al. 2017) and intermediate-old stars (Bagheri et al. 2013;Carrera et al. 2017;Jacyszyn-Dobrzeniecka et al. 2017). Because of their proximity, namely ∼ 50 kpc the LMC (de Grijs et al. 2014) and ∼ 60 kpc the SMC (de Grijs & Bono 2015), they offer a unique opportunity to investigate in great detail the physical effects of a threebody encounter, among the LMC, SMC, and the MW. According to N-body simulations, the expected outcome of these interactions is a wealth of complex and extended tidal features, that may be very difficult to detect because of their very low surface brightness (e.g., Bullock & Johnston 2005;Cooper et al. 2010). Moreover, according to simulations, the origin and the present-day shape of these features require not just tidal effects, but also hydro-dynamical effects of interaction between the MCs (e.g. Wang et al. 2019, and references therein).
In the last decade, several of these structures have been discovered, thanks to modern deep panoramic surveys (e.g, Mackey et al. 2016;Belokurov & Koposov 2016;Pieres et al. 2017;Choi et al. 2018;Mackey et al. 2018;Massana et al. 2020;El Youssoufi et al. 2021). These discoveries aided to decipher the interaction history of the MCs, whose morphology and kinematics have been shaped in the last Gyrs by their mutual attraction and by the gravitational interference induced by the MW (e.g., Besla et al. 2012;. Even more valuable, the Early Data Release 3 (EDR3) of the Gaia mission (Gaia Collaboration et al. 2016, 2021a made it possible to kinematically characterize the already known tidal substructures (e.g., Zivick et al. 2019;Schmidt et al. 2020;Omkumar et al. 2021;Cullinane et al. 2022;James et al. 2021;Piatti 2022), to unveil the response of the MW halo to the LMC (e.g., Garavito-Camargo et al. 2019;Conroy et al. 2021;Vasiliev et al. 2021), and also to discover new stellar streams De Leo et al. 2020;Grady et al. 2021;Petersen et al. 2021).
Amongst the numerous sub-structures revealed in the last years, the arm-like stellar over-density to the North of the LMC, currently known as Northern Tidal Arm (NTA), discovered by Mackey et al. (2016) thanks to the analysis of the public Dark Energy Survey (DES), is particularly prominent. This sub-structure was traced in all its extension by  using red giant branch (RGB) stars selected by means of Gaia Data Release 2 (DR2, Gaia Collaboration et al. 2018). Different authors ran N-body simulations of a MW-LMC interaction to explain the origin of this feature (Mackey et al. 2016;Cullinane et al. 2022). For example,  simulated the evolution of the LMC particles in the last 1 Gyr in the presence of the only MW or SMC, and also under the influence of both galaxies. They also sampled the simulation outcomes with different LMC and SMC mass values. They concluded that the MW is the main actor producing the bending of the Northern part of the LMC disc, but a combination of both the MW and the SMC is needed to originate the ensemble of tidal substructures observed around the LMC. Cullinane et al. (2022) also produced a suite of dynamical simulations of MW-LMC or MW-LMC-SMC interactions, letting them evolve during the last 1 Gyr, and treating each system as a particle sourcing potential in order to consider the response of the MW motion to the LMC. Similarly to , they inferred that the tidal force of the MW is the primary responsible of the LMC Northern bending. It is worth noticing that all these simulations predict a more diffuse twisted stream, located towards the North-East of the LMC, rather than a thin and horizontal arm as observed so-far.
In this paper, we use Gaia EDR3 results to select stars belonging to the intermediate and old population of the LMC to trace the very low-surface brightness regions in the outskirts of the LMC. This procedure allows us to discover a more diffuse stellar structure protruding from the outer disk of the LMC in the North-East direction, right below the NTA and in the same position predicted by existing dynamical models (Mackey et al. 2016;Cullinane et al. 2022).

SAMPLE SELECTION
To carry on a preliminary selection of stars potentially belonging to the LMC, we performed a query to the Gaia database 1 to retrieve astrometric and photometric information from Gaia EDR3 (Gaia Collaboration et al. 2021a). In particular, we confined our research in a sky region defined by the Galactic longitude 200 ≤ l ≤ 360 and Galactic latitude −90 ≤ b ≤ −10, in order to avoid likely MW disc stars which can cause significant contamination to the sample. Moreover, we restricted the query to objects whose PMs are within 3σ from the mean LMC and SMC PMs as measured by Gaia Collaboration et al. (2021b). Following a procedure similar to that described by Gaia Collaboration et al. (2021b), we minimized contamination from MW foreground stars, by requiring that ≤ 3 σ , where and σ are the parallax and the parallax uncertainties, respectively. The latter selection allowed us to exclude sources not compatible with the LMC distance. In addition, we filtered out the stars having errors larger than 0.2 mas yr −1 in each of the two PM components and the sources with RUWE 2 > 1.3 astrometric excess noise sig > 2 to limit our sample to objects with good PMs accuracy and reliable astrometric solutions (see App. A for a reproduction of the query). This query yielded a total of 4, 363, 722 stars.
We used the photometry provided by Gaia to select both RGB and red clump (RC) stars in the colourmagnitude diagram (CMD), using G, G BP , G RP Gaia filters. We focused only on the stars belonging to the intermediate-old populations, because they are known to be distributed at large distances from the LMC centre (e.g., Gaia Collaboration et al. 2021b;El Youssoufi et al. 2021). As a first step, we corrected the photometry for the reddening, to avoid that highly extincted MW disk stars could contaminate our sample. To this end, we adopted the reddening maps by Skowron et al. (2021) that provided extinction values in the inner regions of the LMC and SMC, i.e. ∼ 6 • and ∼ 5 • from the LMC and SMC centre, respectively. Outside these regions, we sampled the reddening maps by Schlegel et al. (1998) by means of the python package DUSTMAPS (Green 2018). We calculated the extinction in the Gaia filters through the relations provided by Wang & Chen (2019). Using the dereddened photometry, we finally selected bona-fide RGB and RC stars by including in our sample all the stars inscribed within the polygon shown in Fig. 1 (see also App. A). This procedure yielded a total of 1, 759, 796 stars, that after a quick visual investigation about their sky positions we deduced are a mix of LMC and SMC stars, with a still not negligible contamination by MW stars (see next section). To further remove the MW noise, it is useful to adopt a new reference frame centred on the LMC, in which the new PMs represent the motion of a star in the plane of the LMC disc. To this aim, we took advantage of the transformation laws reported in Gaia Collaboration et al. (2021b) (see their Eqs. 5 and 6) to define the new coordinates ξ, η, the new PMsξ,η, which represent the coordinates and PMs of a star moving on the flat disc defined by the plane of the LMC 3 . In this new reference system we also calculated the radial and tangential velocities v R , v φ , and the third component of the angular momentum L z . The quantity L z is particularly useful to further clean our sample. In particular, we estimated the average LMC L z,LM C by using all the stars within 8 • from the LMC centre and filtered out all the stars that did not satisfy the condition L z,LM C − 5σ L z,LM C < L z < L z,LM C + 5σ L z,LM C , where σ L z,LM C is the standard deviation of L z,LM C . This cut yielded a total of 1, 452, 155 stars.

SELECTION OF LMC STARS THROUGH MACHINE LEARNING TECHNIQUES
To discern stars belonging to different stellar subcomponents we used a Gaussian Mixture Modelling (GMM) 4 , which estimates the statistics of the Gaussians of each underlying single stellar sub-population and, in turn, evaluates for each star the probability of belonging to these sub-populations. The GMM is an unsupervised clustering technique which needs only one parameter, i.e the number N of finite Gaussian functions that describes the N sub-populations within an overall population. We run the GMM with N = 3 Gaussian to model the three stellar components of our sample: the LMC, the SMC and the MW. As input for the GMM we adopted six parameters 5 : the classical PMs (µ α , µ δ ); the parallax , the radial and tangential velocities defined in the new reference frame v R , v φ ; and the de-reddened colour of the stars (G BP − G RP ) 0 . Finally, we assigned to each stellar sub-component, i.e. LMC, SMC and MW, only stars having at least 99% of probability to belong to one of the three galaxies, while the remaining objects are filtered out. After this step our total sample is reduced to 1, 290, 166 stars.  Table 1 lists the statistics of the three Gaussians derived for each sub-population by the GMM, as well as the number of stars assigned to each component. It can be seen that the MW star residuals were ∼ 1.3%, while the percentage of SMC stars raises to more than 5%, as expected because of the large overlapping of their stars on the CMD. The average PMs of the MCs are in very good agreement with those estimated by Gaia Collaboration et al. (2021b) using a selection of RGB stars, e.g. µ α = 1.77 ± 0.44 mas yr −1 ; µ δ = 0.33 ± 0.63 mas yr −1 for the LMC and µ α = 0.71 ± 0.36 mas yr −1 ; µ δ = −1.22 ± 0.29 mas yr −1 for the SMC, making us confident about the reliability of the GMM technique. Figure 2 shows the allocation of each star to one of the three components, displayed through an orthographic projection centred on the LMC (we used Eq. 1 of Gaia Collaboration et al. 2021b). An inspection of the figure reveals that the GMM technique was able to properly separate the LMC from the MW and the SMC, with two exceptions. The first is the wrong assignment of LMC stars to the MW. These stars are mainly located within the LMC bar and we speculate that the more chaotic kinematics (see e.g. Gaia Collaboration et al. 2021b) of the central bar misleads the GMM so that it incorrectly associates these stars to the MW. Furthermore, the v φ values within 4 • − 5 • to the LMC centre is different with respect to those measured in the outer regions (see e.g., Gaia Collaboration et al. 2021b, their   14). The second is a cloud of stars ( 0.02% of the stars assigned to the LMC), that we call "MW noise" hereafter, located between 20 − 50 • North of the LMC (top-left panel) which we speculate are MW misclassified stars. To confirm this hypothesis, we plotted in Fig. 3 the PMs of LMC stars placed beyond 8 • from the LMC centre (gray points), of the MW noise (green points) and the stars assigned by the GMM to the MW situated in the same position of the MW noise (yellow points). The position and PM distribution of the latter population suggest that the MW noise is likely a tail of the MW sample that is being misclassified because of a significant overlap with the corresponding distributions of genuine LMC stars.
Regarding the stars assigned to the SMC (top right sub-panel in Fig. 2), besides its main body it is also visible the Small Magellanic Cloud Northern Over-Density (SMCNOD) discovered by Pieres et al. (2017), which attests the ability of the GMM at revealing also very low-surface brightness structures (Pieres et al. 2017, estimated µ V = 31.2 mag arcsec −2 for the SMCNOD)

DISCUSSION
In this section we discuss separately the morphology and velocity fields of the outer LMC substructures unveiled by the GMM algorithm.

Morphology
The spatial distribution of the stars selected as LMC members by the GMM is shown in Fig. 4. We only displayed objects located beyond 8 • from the LMC cen- tre to better visualize the outer regions (similar density plots with different choices of the bin size and/or smoothing width are displayed in Fig. 8). The morphology of the LMC towards the Galactic disk (i.e., North-East direction) is strongly disturbed as witnessed by several sub-structures emanating from the disc edge. The most noticeable one is the already introduced NTA which is identified by a red polygon in the figure. It is also possible to discern some elongated features towards the East, up to (X, Y ) (20 • , 0), that we dubbed Eastern Substructure (ES hereafter, highlighted with a magenta rectangle in Fig. 4, which represent the Eastern substructure I (de Vaucouleurs 1955) and the recently unveiled Eastern substructure II (El Youssoufi et al. 2021). In addition, the GMM-selected sample of very likely LMC stars also reveals the ensemble of Southern substructures already discovered by  and Mackey et al. (2018) that we collectively call S-Sub and highlight with an orange polygon. The figure also reveals a diffuse sub-structure protruding from the outer LMC disc and extending at more than 20 deg from the LMC centre and placed in between the NTA and the ES. This feature, indicated with a blue polygon in Fig. 4, represents, as far as we know, a still undiscovered tidal feature, which we call North-Eastern Structure (NES). The NES presents a horizontal overdensity of stars ("Finger" hereafter) at Y 11 • − 12 • which mimics the NTA at about 5 • South of this feature and a more diffuse component which departs from the outer disk at Y 5 • − 7 • deg and connects itself with the Finger at about (X, Y ) (15 • − 20 • , 10 • ). Figure 5 shows the CMD of the stars belonging to the NTA, the NES, the ES and the MW noise (highlighted with different colours), superposed with the CMD of LMC stars at its North-East side and beyond 9 • . The NTA (top-left panel) and the NES (top-right panel) have an RC-RGB similar to the underlying population of the outer LMC disc. This occurrence strengthens the hypothesis that the NES is composed by stars tidally stripped from the outer LMC disc, as the NTA. To be more quantitative, the ratio (f ) between the number of RC (indicated with the black polygon in the figure) and RGB stars for the LMC outer disc, the NES and the NTA are very similar, being f = 0.34 ± 0.01, f = 0.29 ± 0.05 and f = 0.28 ± 0.05, respectively. The ES, instead, shows an almost total lack of RC stars.
Such lower fraction for the ES might be due to issues with the completeness and/or reddening estimates, however investigating this occurrence is beyond the purpose of this work. Finally, the MW noise (bottom right panel) contains stars mostly having G 0 > 18, reinforcing the hypothesis that they are MW misclassified stars.
To further strengthen the reality of the NES, we point out that some hints of its presence were already visible in other works that exploited the Gaia catalogue. For example, the Finger is barely visible in the density maps constructed with Gaia EDR3 in Gaia Collaboration et al. (2021b) (their fig. 17), but the NES is more visible in fig. 2 of Grady et al. (2021), which shows a selection of 226,119 RGB stars in both the LMC and SMC from the Gaia DR2. The main scope of their work was to obtain photometric metallicity estimations and an overall picture of the PMs of their sample of stars rather than discover new undetected substructures in the faint outer regions of the LMC, therefore they also did not need to remove the minor MW component from their sample. The GMM adopted in this work cleaned the RGB stars from the MW contamination, yielding an almost pure catalogue of LMC stars, hence particularly enhancing the NES visibility.
To investigate the origin of the NES we inspected the recent results obtained using N-body simulations of the interaction between the LMC, SMC and MW. Mackey et al. (2016) ran N-body MW-LMC simulations to understand the physical origin of the NTA. They demonstrated that the tidal interaction with the MW can warp and distort the outer LMC disk, generating two spirallike patterns, one twisted towards the North-East and the other towards the South-West direction. Their results strongly indicate that the NTA originated by the tidal interaction of the LMC with the MW. Similar structures have been predicted by  N-body models. In particular, simulations of the LMC-SMC-MW or of the MW-LMC interactions, under some initial conditions 6 foresee both a leading and trailing spirals, stretched in the Southwest-Northeast direction (see their Fig. 4). They speculated that the NTA originated from the gravitational interaction with the MW, but a combined influence of the MW and SMC is necessary to explain all the features observed in the South of the LMC. Interestingly, these models also forecast, as a distinct and separated feature from the NTA, the bending of the entire Northern outer disc of the LMC caused by the MW gravitational pulling. Therefore, besides a thin horizontal stellar over-density, one should observe a more diffuse and wide stellar structure which should represent the distorted Northern LMC disk. However, observational works until now did not disclose any diffuse stellar over-density in the North, but only the thin NTA. Morphologically, the NES is in the position expected from the quoted dynamical models and we can hypothesize that it might have had the same origin of the NTA. In this scenario it would represent the tidal distortion induced by the MW to the outer North-East of the LMC disk.
To be complete it is worthwhile to underline that some models of the LMC interacting with the SMC and the MW do not expect the presence of a diffuse bended outer structure. For example, Besla et al. (2016) brought to light the existence of stellar arcs and multiple stellar arms to the North of the LMC. From hydrodynamical simulations of the MCs either when evolved in isolation or under the gravitational influence of the MW, they concluded that repeated close interactions with the SMC are the main cause of the many stellar substructures observed in the Northern periphery, included the NTA. However, their models were based on the observations that a Southern counterpart of the features revealed to the North does not exist. These Southern substructures (S-SUB in Fig. 4) were disclosed afterwards by Mackey et al. (2018) and . Therefore, because of the asymmetric features revealed only to the North, Besla et al. (2016) disfavor a formation scenario driven by global Milky Way tides. This is at odds with the more recent models carried out by  and Cullinane et al. (2022) which are based on more recent and complete data on the features unveiled around the LMC.

Velocity field
If the NES and NTA had the same origin, they should also show a coherent velocity field. In the context of the MagEs survey, Cullinane et al. (2022) analyzed the kinematics of the NTA through spectroscopic measurements of different fields centred along the arm, in conjunction with Gaia EDR3 PMs. From a combined study of NTA metallicity and overall kinematics, they definitely confirmed that the NTA is made of LMC outer disc material tidally distorted by an external gravitational potential. They also carried out a suite of N-body simulations, ob- taining that the LMC embedded in a heavy-MW potential (mass MW ∼ 10 12 ) represents the best qualitative match between models and observations. However, their best model shows a diffuse stellar component connecting the LMC disc with the NTA, which resembles the NES. From the dynamical point of view, the Cullinane et al. (2022) simulations are not able to reproduce the measured kinematics of the NTA. Indeed, while the measured tangential velocity is quite well replicated by the best model, the radial and vertical velocities are inconsistent with stars in an equilibrium disc and very different from what is expected by a solely LMC-MW interaction.
To further investigate this point, Fig. 6 (top row) shows the LMC v φ (left panel) and v R (right panel), calculated as median of the velocity of the stars within pixels of 0.7 • × 0.7 • , and their dispersions (bottom row). Both these diagrams seem to suggest a similar velocity pattern for the NTA and the NES, constituted by strong velocity gradients, with their velocity decreasing either going from the ES towards the NTA or moving along the NTA towards the West direction. According to the velocity dispersion maps in the lower panels of Fig. 6 (the one in σ R , in particular), NES and NTA appear dynamically colder than the adjacent edge of the LMC disc (around (X, Y ) (5 • , 12 • ) or than other tidal structures, such as S-SUB, which is very hot (i.e. σ R ≥ 70 km s −1 ). However, some of the pixels beyond Bottom: Velocity dispersion of v φ (left panel) and vR (right panel). Each pixel is 0.7 • × 0.7 • in size. The coloured polygons are the same as in the right panel of Fig. 4 10 • − 15 • include few stars, causing the outer velocity map of the LMC to be noisier and more scattered. It is also visible from the v R plot that the NTA and the Finger present negative radial velocities (top right panel) of the order of −50 km s −1 , similar to those estimated by Cullinane et al. (2022), with peaks of −80 km s −1 .
To explain the origin of the peculiar observed radial velocities of NTA's stars Cullinane et al. (2022) speculated that they could have been caused by recurring close passages of the SMC, happened more than 1 Gyr ago, beyond the time span of their models, which cover only the last Gyr. According to their simulations, any close passage of the SMC in this period was not able to alter significantly the kinematics of the NTA. In this context, the Finger might have been affected by past nearby passages of the SMC as well, hence it might be connected to the NTA and physically distinct from the rest of the NES. However, the velocity maps of the NTA and NES displayed in Fig. 6 do not allow us to establish whether the NES-NTA is a single highly disturbed sub-structure or is made by two uncorrelated physical sub-structures originated at different times. An additional characteristic visible in Fig. 6 is the large positive radial velocity shown by ES stars, with peaks of v R 80 km s −1 , possibly indicating that the ES is in an ongoing stripping process.
An overall picture of the bulk motion of the outer regions, corrected for the solar reflex motion, is depicted in Fig. 7. At first glance the plot suggests that the NTA is being pulled back by the gravitational potential of the LMC towards equilibrium, as already visible in Fig. 6. The overall kinematics of the ES, and more generally of the outer East-North-East side of the LMC, indicates that this side is strongly moving away from the LMC. The NES, instead, shows a more chaotic kinematics, as some regions have a motion more similar to the ES while others have an opposite velocity direction. In the same figure the highly disturbed kinematics in the S-SUB is also visible. The high perturbations induced by the MW is testified by the strong deviations from the rotational motion of the East LMC side with respect to its West side. However, we caution the reader that the LMC inclination may affect the velocities discussed above and the consequent results. A quantitative discussion of this scenario could be provided only when the velocity field is drawn also with spectroscopically inferred radial velocities. Recently, also Choi et al. (2022) revealed large PM residuals in the Southern disc with a sample of RC and RGB stars collected from Gaia EDR3, even though their analysis was limited to the inner 6 • from the LMC centre. They compared the observed LMC PMs with the outcome of numerical simulations of a LMC in interaction with the MW and SMC, putting constraints on the lookback time the last close encounter between the MCs happened (i.e., t < 250 Myrs ago) and on its impact parameter (i.e., ≤ 10 kpc.) 7 .
An additional interesting feature in the North-East side of the LMC was recently discovered by Petersen 7 These constraints are very similar to those inferred by Zivick et al. (2018) et al. (2021) using the RR Lyrae variables as kinematics tracers. Indeed, they detected a tidally stripped stream, expected to be the outcome of the MW tidal force, in the same direction of the features we have discussed in the North-East side of the LMC. Petersen et al. (2021) probed only regions beyond 20 − 25 • from the LMC centre, and revealed the stream signal up to ∼ 30 • (see their fig. 1), finding also likely stream members in the Northern hemisphere at about 70 • from the LMC centre. In this context, we speculate that the ensemble of stellar streams in the North-East outer disc could represent the inner part of the outer stream disclosed by Petersen et al. (2021), which extends up to the Galactic disc.
Dynamical simulations of LMC-SMC-MW interactions, going back to at least 2-3 Gyr, are essential to explain the peculiar kinematic features observed in the North-East LMC, such as the strong velocity gradients along the NES-NTA sub-structures.

SUMMARY
In this work we exploited the exquisite multidimensional data provided by Gaia EDR3 to investigate the outer low surface brightness substructures of the LMC. We selected LMC RGB and RC stars based on their position on the CMD and filtered out all the stars having angular momentum beyond 5 σ that of the average LMC disc. We run the unsupervised clustering technique GMM with 6 input parameters to separate the remaining stars into the three main stellar components, the MW, LMC and SMC, and we limited our analysis to highly-likely LMC members (P > 99%).
The GMM was able to recover the outer sub-structures discovered in the last years in the LMC, such as the NTA, the S-SUB or the ES. In addition, this procedure allowed us to detect an unknown diffuse tidal structure at its North-East side extending up to ∼ 20 • from the LMC centre, the North-East Structure (NES). The NES fills the gap between the outer LMC disk and the other known structures in the North-East of the LMC, namely the ES and the NTA. The presence of tidally distorted LMC stars at the location of the NES was expected based on N-body simulations of a MW-LMC interaction Cullinane et al. 2022) as a consequence of the tidal stress induced by the MW on the LMC. This occurrence might induce us to think that the NES has the same origin as the NTA. From the kinematic analysis based on the Gaia EDR3 PMs, these tidal stellar components reveal a strong velocity gradient towards the outer LMC, but the present dataset is not sufficient to assess whether or not the NES is physically correlated to the NTA. Overall, the North-East side of the LMC displays radial velocities not consis-tent with a disc in equilibrium, showing large negative velocities for the NTA-Finger and positive for the ES. This occurrence indicates that the outer LMC disc has been strongly disturbed in the past few Gyrs. Spectroscopic follow-up to measure accurate radial velocities and more detailed and extended in time N-body simulations are requested to assess the origin of the NES and of the other sub-structures in the North-East side of the LMC.
We warmly thank the anonymous Referee for the comments that helped us to improve the manuscript. This work has been partially supported by INAF through the "Main Stream SSH program" (1.05.01.86.28). M.G. and V.R. acknowledge support from the INAF fund "Funzionamento VST" (1.05.03.02.04).

B. DENSITY MAPS OF THE OUTER LMC
In this section we show the density map of the LMC as depicted in Fig.4, with different choices of bin size and width of the kernel function. To enhance the small substructures we adopted narrower kernel function, while to have a global view of the outer LMC we set a larger bin size and shallower kernel function.