Offset of M54 from the Sagittarius Dwarf Spheroidal Galaxy
Zhaozhou An Matthew G. Walker Andrew B. Pace


firstpage–lastpage
[
    [
    25 April 2024
=================




    We present results from simultaneous modeling of 2D (projected along the line of sight) position, proper motion and line-of-sight velocity for Gaia- and APOGEE-observed stars near the centre of the Sagittarius (Sgr) dwarf spheroidal galaxy.  We use a mixture model that allows for independent sub-populations contributed by the Sgr galaxy, its nuclear star cluster M54, and the Milky Way foreground.
    We find an offset of 0.295± 0.029 degrees between the inferred centroids of Sgr and M54, corresponding to a (projected) physical separation of 0.135± 0.013 kpc.  The detected offset might plausibly be driven by unmodelled asymmetry in Sgr's stellar configuration; however, standard criteria for model selection favour our symmetric model over an alternative that allows for bilateral asymmetry.  We infer an offset between the proper motion centres of Sgr and M54 of [Δμ_αcosδ,Δμ_δ]=[4.9, -19.7] ± [6.8, 6.2] μas yr^-1, with magnitude similar to the covariance expected due to spatially-correlated systematic error.  We infer an offset of 4.1± 1.2 km s^-1 in line-of-sight velocity.  
    Using inferred values for the systemic positions and motions of Sgr and M54 as initial conditions, we calculate the recent orbital history of a simplified Sgr/M54 system, which we demonstrate to be sensitive to any line-of-sight distance offset between M54 and Sgr, and to the distribution of dark matter within Sgr.  



    galaxies: dwarf – galaxies: individual: Sgr dSph – globular clusters: individual: M54 – Galaxy: structure






§ INTRODUCTION



Observations reveal that a significant fraction of galaxies host dense nuclear star clusters (NSCs) at their centres <cit.>.  However, the
formation of NSCs is still a debated question.  Possible explanations range from in situ star formation at the galactic centre
<cit.> to the infall of globular clusters (GCs), whose orbits decay and spiral to the centre due to
dynamical friction <cit.>.

Some NSCs are found to host metal-poor stars <cit.>, which can be naturally explained by
old, metal-poor GCs bringing those stars into the nucleus.  Deficits of bright non-nuclear GCs near the central regions of elliptical galaxies also seem to indicate
that the central nuclei formed via the orbital decay of massive clusters
<cit.>. In early-type galaxies with stellar mass ≲ 10^9 M_⊙, the NSCs and GC systems have similar occupation distributions and comparable total masses <cit.>, consistent with NSC formation by infalling GCs.  

However, many NSCs also contain young stellar
populations <cit.>, suggestive of in situ star formation within the host galaxy.  Thus, it seems that both the
GC infall and in situ star formation contribute to NSC formation.  Indeed some NSCs display both young and old stellar populations with different morphology
and kinematics, which indicates  different origins for different
components and multiple star formation episodes <cit.>. For example, <cit.> show
that it is necessary to consider both stellar dynamical mergers and in situ
star formation to reproduce observations of the NSCs at the
centres of NGC 422 and M33. <cit.> use a semi-analytical
galaxy formation model that includes both dynamical-friction-driven migration
of stellar clusters and star formation triggered by infalling gas to simulate
the co-evolution of NSCs and central black holes, finding that in situ star formation contributes a
significant fraction (up to ∼ 80%) of the total mass of NSCs.


Considering the various pathways of NSC formation, it is useful to make 
 detailed comparisons between the spatial and chemodynamical properties of NSCs and their host galaxies.  Notably, some works find evidence for offsets between the centres of some NSCs and the centres of their hosts.   
For example, <cit.> search for off-centre nuclei in 78
nucleated dwarf elliptical galaxies of the Virgo Cluster, finding that roughly
20% of the sample are significantly off-centre.
<cit.> study 100
early-type galaxies of the Virgo Cluster, identifying five candidates with nuclei offset from the host galaxy photocentre; they note, however, it is possible that some of 
these cases are due to non-nuclear globular clusters that happen to be projected close to the galaxy photocentre.  Such offsets, if real, may indicate oscillations of the NSC that may help
constrain, e.g., the mass of a central black hole <cit.> and/or other physical processes like counterstreaming instability <cit.>.

The nearest NSC, Messier 54 (M54 also known as NGC 6715) within the Sagittarius (Sgr) dwarf spheroidal (dSph) galaxy, offers a unique opportunity to observe detailed spatial and chemodynamic properties of individual stars within a NSC and its local environment.  Originally discovered in a line-of-sight velocity survey of the Galactic bulge <cit.>, Sgr lies at a distance of ∼ 26 kpc from the Sun,  ∼ 18 kpc from the Galactic centre, and is the progenitor of the prominent Sagittarius stream, stellar debris from the ongoing process of Sgr's tidal disruption within the Milky Way <cit.>.  The central region of Sgr, which includes M54, displays a complex formation history, with various photometric and spectroscopic studies identifying old (≳ 10 Gyr), intermediate-age (∼ 4-6 Gyr) and young (≲ 3 Gyr) stellar populations with a corresponding range of metallicity from -1.8≲ [Fe/H]≲ +0.6 <cit.>.  

This complexity extends to the observed stellar kinematics.  The old, metal-poor stars have a line-of-sight velocity dispersion profile that declines steeply with cluster-centric radius <cit.>, then gradually rises near M54's nominal tidal radius <cit.>.  Over the same region, the relatively metal-rich stars display approximately constant velocity dispersion, as is characteristic of the Milky Way's dwarf spheroidal satellite galaxies <cit.>.  Even outside the nucleus, the main body of Sgr itself exhibits complex internal stellar kinematics, with a central bar-like structure that connects to `tails' of escaping stars, a bound, rotating inner core,  and apparent expansion along the long axis that characterises an overall triaxial morphology   
<cit.>.  

Given this complexity, it is difficult to disentangle M54 from its Sgr host.  Analysing the spatial distribution of a large photometrically-selected sample of M giant stars from the Two Micron All Sky Survey,  <cit.> find an overdense `cusp' at the position of M54; however, they argue that this feature belongs to Sgr and not to M54, whose metal-poor stars are too blue to be included in the M giant sample.  Analysing metal-rich and metal-poor populations selected from optical photometry, <cit.> conclude similarly, finding that the metal-rich stars (presumed to belong to Sgr) exhibit an overdense cusp at the same location as the cusp seen in metal-poor stars (presumed to belong to M54).  However, a MUSE spectroscopic survey of M54 by <cit.> finds that the metal-rich stars can be separated into young (∼ 2 Gyr) and intermediate-age (∼ 4 Gyr) populations, associating the more compact, younger subpopulation with late star formation within M54, and the more extended, intermediate-age sub-population with the Sgr host.

If the old, metal-poor component of M54 is an example of a GC that fell to the centre of its host galaxy, then its properties may also help to constrain the spatial distribution of dark matter within Sgr.  Various dynamical and N-body studies demonstrate that the efficiency and outcome of dynamical friction depend on the spatial distribution of dark matter within the host galaxy <cit.>.  If the host dark matter halo has the central `cusp' that characterises cold dark matter halos <cit.>, dynamical friction operates efficiently to drag a massive GC to the centre in a fraction of a Hubble time; if instead the host halo has a central `core' of uniform density, the infalling cluster tends to stall near the core radius, where the halo becomes effectively `bouyant' <cit.>.  Using N-body models to simulate specifically the infall of M54 within Sgr, <cit.> find that for a variety of initial conditions, if the Sgr host has an NFW-like cusp, then M54 sinks efficiently (in ≲ 3 Gyr) to the very centre of Sgr, where it remains virtually at rest, with a systemic velocity offset of ≲ 1-2 km s^-1.  These results are broadly consistent with those of <cit.>, who also find that a more (than NFW) steeply cusped dark matter halo in Sgr would result in a NSC that (internally) rotates faster and is morphologically flatter than M54. 


Thus the observed properties of the Sgr/M54 system impact our understanding not only of NSC formation, but also the nature of dark matter.  In any case, the degree to which M54 and the centre of Sgr are either distinct or coincident in phase space remains an interesting open question.  The answer is complicated for several reasons.  First, as mentioned above, the central regions of Sgr display a mixture of several stellar populations with different ages, metallicities, structural parameters and kinematics, and the assignment of a given population to one or the other object can be ambiguous.  Second, published measurements of centroids, distances and systemic velocities often pertain to just one of Sgr or M54.  For example,  the infrared M-giant sample of <cit.> is relatively insensitive to the metal-poor population in M54; conversely, Hubble Space Telescope and the MUSE observations of M54 typically lack the sky coverage to map Sgr 
 over a significant fraction of its core region  <cit.>.  As a result, comparisons of published measurements for Sgr and M54 must contend with systematic errors (e.g., zero-point offsets) introduced by the use of different data sets.  

A notable exception is the recent study by <cit.>, which uses the Gaia (DR2) catalogue, supplemented with training data from spectroscopic (line-of-sight velocity) catalogues, to construct a catalog of 6D phase-space coordinates for a final sample of ∼ 1.2× 10^5 stars within ∼ 6^∘ of the Sgr centre.  With this catalogue, they iteratively measure the centre-of-mass coordinates of a Sgr sample that excludes all stars within the fiducial tidal radius of M54, and separately the center-of-mass coordinates for the M54 sample that includes only the stars excluded from the Sgr sample. 
 They find a statistically significant offset only among the plane-of-sky centroids, which they discount as likely driven by systematic error due to  incompleteness that  results from the Gaia scan pattern that is apparent in their sample (their Figure 4).  

Here we take an approach that is complementary to that of <cit.>.  Rather than adopting hard cuts in colour/magnitude or sky position to separate M54 and Sgr samples, we use the mixture models to analyse the 5D distribution of sky position, proper motion and line-of-sight velocity, the phase-space coordinates for which homogeneous measurements from Gaia EDR3 <cit.> and/or the Apache Point Observatory Galactic Evolution Experiment <cit.> are available for large numbers of stars in the Sgr.  The mixture modelling lets us obtain simultaneous estimates for the 5D centres of the Sgr and M54 populations, providing a direct inference about any offset in these dimensions of phase space.  Any detected offset can then be used to inform subsequent models for the formation and evolution of M54 as the NSC within Sgr.  Furthermore, since our mixture models do not operate separately on pre-selected M54 and Sgr samples, the result can be used to infer, rather than assume, the colour/magnitude distributions traced separately by the two objects.  

We describe the data selection procedure in Section <ref>. In Section <ref>, we describe the construction of the mixture model.  In Section <ref>, we present results from our modelling. In Section <ref>, we compare our results to previous work and discuss the robustness of our results, the colour/magnitude distribution of the M54, and the recent orbital history of the Sgr/M54 system.  In Section <ref>, we summarise our conclusions.




§ DATA



This work is based on data presented in the Early Data Release 3 <cit.> of the Gaia mission <cit.>, and on line-of-sight velocity data presented by APOGEE Data Release 17 <cit.>.  
Gaia EDR3 consists of astrometry and photometry for 1.8 billion sources brighter than magnitude G≈ 21.  Compared with Gaia DR2, Gaia EDR3 improves the parallax precision by 30% and increases proper motion precision by a factor of 2. At magnitude G=17, it provides typical uncertainties of 0.07 mas yr^-1 in proper motion, 0.07 mas in parallax angle, 0.05 mas in position, 1 mmag in G-band magnitude, 12 mmag in G_BP magnitude, and 6 mmag in G_RP magnitude <cit.>. 
APOGEE DR17 is the continuation of the Apache Point Observatory Galactic Evolution Experiment <cit.> and Data Release 17 of APOGEE-2/SDSS-IV  <cit.>.  APOGEE spectroscopic catalogues include stellar-atmospheric parameters and line-of-sight velocities measured from twin, multi-plexed, near-infrared, high-resolution spectra covering both the northern and southern sky <cit.>. 
 The extracted line-of-sight velocities have a median reported uncertainty of 0.015  km s^-1. 




 §.§ Gaia EDR3 Data Selection



We begin by considering all sources from Gaia EDR3 that lie within 4^∘ of a fiducial center of M54 <cit.>.  Within this region, one immediately obvious systematic effect is the tendency for sources near the centre of M54 to have relatively poor astrometric solutions due to crowding.  Figure <ref> demonstrates that application of the standard filter for selecting astrometrically well-behaved sources, <2 <cit.> removes almost all sources within ∼ 1 arcmin, or ∼ 2 half-mass radii of M54 (r_h=0.43± 0.08 arcmin; ).  To avoid this 
selection bias, we mask the entire region within ∼ 0.03^∘ (1.8 arcmin) of the fiducial centre of M54. We note that special Service-Interface Function (SIF) images of dense area, which includes the centre of the Sgr dSph, are acquired in the sky-mapper CCDs <cit.>, and that can be used to analyse the centre of the Sgr dSph. One example is the data for ω Centauri in <cit.>. However, these special SIF data are not part of Gaia DR3.

Next, following <cit.>, we impose a magnitude cut at G<17.3 in order to ensure high-quality astrometry and homogeneous coverage over the region of interest.  Then we cut on proper motion based on a data-driven filter as shown in Fig <ref>.
We make additional cuts on parallax to remove stars with distance less than 20 kpc and apply data quality flags.  Our complete initial selection uses the following criteria:



0.03 ≤√(x^2+y^2)≤ 4

<17.3,

<2,

parallax<0.05+3 × parallax_error,

pmra∈ [-3.12, -2.22], pmdec∈ [-1.91, -0.91]

ipd_frac_multi_peak < 2

ruwe < 1.3

duplicated_source =  False

| C^⋆|<3σ_C^⋆.



Here, (x,y) are rectilinear coordinates obtained via gnomonic projection of the Gaia-measured sky coordinates, assuming an origin at the previously-published center of M54 <cit.>.  C^⋆ is the corrected BP and RP flux excess factor, defined in Equation (6) of <cit.> and σ_C^⋆ is the
1σ scatter, calculated according to Equation (18) of
<cit.>.  Using the `astrometric fidelity' parameter of <cit.>, all of the sources in our final sample have fidelity>0.5 <cit.>, indicating good data quality.  The Gaia EDR3 photometry is corrected for extinction using the code from Appendix A of <cit.>. The extinction is calculated using Equation 1 and Table 1 from <cit.>.










 §.§ APOGEE DR17 Data Selection

From APOGEE DR17 we select stars within a 4-degree circle centred at (α,δ)=(283.76353, -30.477006).  We remove stars for which any of the following conditions are met :  (bit position 23 in ) is set to 1,  is nonzero, or  (bit position 4 in ) is set to 1. The APOGEE DR17 catalogue is cross-matched with Gaia EDR3, allowing us to apply to the APOGEE sample the same proper motion selection that we adopt for the Gaia EDR3 sample (Sec <ref>).  Our selections reduce the number of stars in our APOGEE DR17 sample from 6890 to 683.
Fig. <ref> shows the spatial distribution of stars in the selected APOGEE sample.  Since APOGEE DR17 has only partial coverage in this sky area, we cannot directly combine the APOGEE and Gaia samples.  Therefore, below, we choose to analyse the line-of-sight velocity data separately from the sky position and proper motion data.  
(Sec <ref>).







 §.§ CMD mask

The samples selected according to the above criteria will have contamination from Milky Way stars that belong to neither M54 nor Sgr.  In order to minimise this contamination, we construct an additional mask in the colour/magnitude diagram (CMD), combining a
mask derived from the data and a mask derived from isochrone templates.

To construct the CMD mask from the data, we consider two annuli that are centred on M54, the first with inner and outer radii of 0.03, 0.2, respectively, and the second
with radii of 2.0, 4.0.  Then, from the selected sources within each annulus, we calculate  histograms in colour-magnitude space, with 40 uniform bins from -0.2 to 2 in the BP-RP colour
dimension and 40 uniform bins from 13 to 19 in the G magnitude dimension (only for the purpose of constructing the CMD masks, we relax the magnitude limit of G< 17.3 that we apply when analysing the samples).  Both 2-D histograms are divided by the annulus area, so that the result indicates the density of stars.  We then subtract the histogram corresponding to the outer annulus from the histogram corresponding to the inner annulus.  Since the outer annulus should contain a larger fraction of Milky Way contaminants, we expect this difference to be more positive in colour/magnitude bins that are less contaminated by the Milky Way.  We mask sources within colour/magnitude bins for which the difference histogram has value <20 stars per degree2. In order to decrease the effect of noise, we perform a `closing' operation, using a  2 × 2 matrix of ones as the structuring element on the CMD mask. 

To construct the CMD mask from isochrone templates, we consider three sets of theoretical isochrones that are motivated by the age/metallicity relations derived from MUSE spectroscopy of M54 by <cit.>: young metal-rich, with [Fe/H]=-0.04 and age between 1.4-2.8 Gyr, 2) intermediate-age metal-rich, with [Fe/H]=-0.29 and age between 1-7.8 Gyr, and old metal-poor, with [Fe/H]=-1.41 and age between 9.4-15 Gyr.  For each set, we obtain the corresponding isochrones in Gaia passbands from the Dartmouth Stellar Evolution Database (DSED) <cit.>.  We sample the age ranges with a step size of 0.2 Gyr.  Shifting each isochrone for a fiducial distance modulus of 17.27 <cit.>,  we then mask all CMD bins (using the same binning scheme adopted above for the data-derived mask) that do not contain any of the sampled isochrones. 
 

Finally, we combine data-derived and isochrone-derived masks by masking only those bins that are masked in both schemes. 
 Fig <ref> shows the data-derived, isochrone-derived and combined CMD masks.  Application of the combined mask reduces the number of stars in our selected Gaia EDR3 sample from 44337 to 34265. 
 





§ MIXTURE MODEL FOR THE SGR/M54 SYSTEM



We assume the sample of unmasked stars includes contributions from Sgr, M54 and the Milky Way.  First we construct and fit a mixture model for the joint 4D distribution of Gaia-observed positions and proper motions.




 §.§ Projected position



For M54 and Sgr, we assume that the 2D positions, 

    s=[ x; y ],

of stars are drawn independently from analytic surface density profiles that allow for elliptical symmetry about centers specified by

    s_0=[ x_0; y_0 ].

For ellipticity ϵ≡ 1-b/a specified by semi-major and semi-minor axes a and b, and with semi-major axis pointing along position angle θ_0 (increasing northward of the zero-point at due East), a star's `elliptical radius', R_e, relative to the center satisfies 

    R_e^2 = R^⊤R=( s-s_0)^⊤R_θ^⊤[         1         0;         0 1/(1-ϵ)^2;           ]R_θ( s-s_0),

where  

    R_θ=[  cos(θ)  sin(θ); -sin(θ)  cos(θ) ]

is the rotation matrix.  

We assume that within the observed region, the 2D spatial distribution of stars belonging to M54 is drawn from the probability distribution specified by <cit.>, 

    f_s,King(s|Θ_King )=k_King(1/√(1+(R_e/r_c)^2)-1/√(1+(r_t/r_c)^2))^2,

with Θ_King≡{x_0,y_0,r_c,r_t,ϵ,θ} the vector of free parameters that specify the centroid, `core' and `tidal' radii, ellipticity and position angle, respectively.  The normalisation constant k_King guarantees that ∫_ field f_s(x, y) dx dy = 1, where the integration integrates only over the sampled field.

 We assume that the 2D spatial distribution of stars belonging to Sgr is drawn from the probability distribution specified by <cit.> 

    f_s,Plummer(s|Θ_Plummer)=k_Plummer/(1+R_e^2/r_p^2)^2,

with Θ_Plummer≡{x_0,y_0,r_p,ϵ,θ} the vector of free parameters that specify the centroid, Plummer scale radius, ellipticity and position angle, respectively, and k_Plummer the normalisation constant.

We assume the Milky Way contaminants are drawn from a probability distribution specified by a first order polynomial function, 

    f_s,n(s|Θ_L)=k_L|1 +  b_0 x + b_1 y |

with Θ_L≡(b_0, b_1) the vector of free parameters that specify the polynomial coefficients, and k_L the normalisation constant.





 §.§ Proper motion

For each of the M54, Sgr and Milky Way components, we assume that the proper motions, 

    μ=[ μ_x; μ_y ],

are drawn independently from bivariate normal probability distributions with centers

    μ_0(x,y |A)=[ μ_0_x; μ_0_y ]+  A[ x; y ],

where 

    A=[ a_1 a_2; a_3 a_4 ]

specifies dependence on projected sky position, and covariance matrix

    Σ=R_ψ^⊤SR_ψ,

where

    S≡[ d_0^2     0;     0 d_1^2 ]

and 

    R_ψ≡[  cos(ψ) -sin(ψ);  sin(ψ)  cos(ψ) ],

with ψ corresponding to the position angle of the long axis of the proper motion ellipsoid (increasing from the positive direction of proper motion in right ascension direction toward  the negative direction of proper motion in declination direction).

For a given component, then, the observed proper motion has probability density

    f_μ(μ|Θ_μ)=𝒩(μ|μ_0,Σ+Σ_0)/C,

where Θ_μ≡{μ_0_x,μ_0_y,a_1,a_2,a_3,a_4,d_0,d_1,ψ} specifies free parameters and Σ_0 is the covariance matrix of the observed proper motion coordinates  (adopted from the Gaia EDR3 catalog), and 

    C≡∫_μ_x_1^μ_x_2∫_μ_y_1^μ_y_2𝒩(μ|μ_0,Σ+Σ_0) dμ_x dμ_y

is a normalizing factor that accounts for the finite proper motion selection window (Section <ref>), chosen to ensure ∬ f_μ(μ|Θ_μ)  dμ=1.  



 §.§ Line-of-sight Velocity



We assume that the line-of-sight velocities, v, sampled by APOGEE are drawn independently from a mixture that includes contributions from Sgr, M54 and the Milky Way, all of which are assumed to have velocities that follow univariate normal probability distributions.  That is, for each of the three components we assume the APOGEE-observed velocities are drawn from 

    f_v(v |Θ_v) = 𝒩(v | v_0 + e_x x + e_y y,σ^2+δ_v^2),

where the free parameters in Θ_v specify a mean velocity, v_0+e_xx+e_yy, that can vary linearly with sky position and has central value v_0, and intrinsic velocity dispersion σ.   δ_v is the observational error associated with the velocity measurement (adopted from the APOGEE catalogue).  In practice, we assume e_x=e_y=0 for the M54 and MW components, allowing for a non-zero line-of-sight velocity gradient only in Sgr.




 §.§ Likelihood Functions


With the normalized probability distributions specified for observed position and proper motion coordinates, the 4D observation of position and proper motion coordinates has joint probability density

    f_s,μ(s,μ|Θ)= ∑_i=1^N_ compϕ_i f_s(s|Θ_s_i) f_μ(μ|Θ_μ_i),

where Θ specifies all free parameters and the sum is over the N_ comp different mixture components, each of which has contribution weighted by ϕ_i, the fraction of the sample that is contributed by the i^th component.  For each of the N_ comp=3 components, we specify ϕ_i in terms of parameters p_i, with the latter defined according to 



    ϕ_i = 
                p_i (1-∑_j=1^i-1ϕ_j )    = p_i∏_j=1^i-1(1-p_j)     i < N_ comp
    
                1-∑_j=1^i-1ϕ_j    = ∏_j=1^i-1(1-p_j)     i=N_ comp


As described in Section <ref>, for different components we adopt different analytic functions for the probability density of stellar positions.  For M54 we use f_s(s|Θ_s)=f_s,King(s|Θ_King); for Sgr we use f_s(s|Θ_s)=f_s,Plummer(s|Θ_Plummer); for the Milky Way we use f_s(s|Θ_s)=f_s,n(s|Θ_L).  

For all three components, we assume that the observed proper motions follow the general form of f_μ(μ|Θ_μ) given by Equation <ref>.  For both M54 and the Milky Way, we assume that the proper motion centre, μ_0, is independent of the spatial position within the sampled region, with a_1=a_2=a_3=a_4=0.  For Sgr, following <cit.>, we allow for a linear dependence of the mean proper motion on projected sky position, letting the components of A be nonzero.

Assuming the observations of different sources are uncorrelated, the 4D Gaia sample of N=34265 stars has likelihood

    ℒ_4=∏_k=1^N f(s_k,μ_k |Θ)


Similarly, the 1D Apogee sample of N_v=683 stars has likelihood

    ℒ_1 = ∏_k=1^N_v∑_i=1^N_ comp q_k,if_v(v_k |Θ_v,  i),

where the weight

    q_k,i = c_i f_s,μ(s_k,μ_k|Θ_i) /∑_j=1^N_ comp c_j f_s,μ(s_k,μ_k|Θ_j)

is the probability of the k^ th star's membership in the i^ th component, given the Gaia-observed position and proper motion parameters Θ_i of the 4D model for the i^ th component.  The three parameters c_1, c_2 and c_3 represent the mixing fractions within the spectroscopic sample. 
 We use the analogue of Equation <ref> to specify the c_i in terms of two independent free parameters m_1, m_2, requiring that ∑_i=1^N_ compc_i=1.




 §.§ Inference

From the 4D and 1D data sets, denoted D_4 and D_1, respectively, we are interested in inferring posterior probability distributions

    f_4(Θ| D_4)=ℒ_4 π(Θ)/∫ℒ_4 π(Θ) dΘ
 
and 

    f_1(Θ_v | D_1)=ℒ_1 π(Θ_v)/∫ℒ_1 π(Θ_v) dΘ_v,

where π(Θ) and π(Θ_v) are priors and the denominator on the right-hand side of each equation is the marginalized likelihood, or `evidence'.  In order to estimate f_1(Θ_v| D_1), when evaluating the component membership probabilities in Equation <ref> we hold 4D model parameters Θ fixed at the maximum a posteriori (MAP) value obtained when estimating f_4(Θ| D_4). 
 Table <ref> lists the priors that we adopt for all model parameters.  For all parameters except M54's core radius and those that specify member fractions, we adopt uniform priors over finite ranges.  For the core radius of M54, since our sample excludes the core region (Figure <ref>), we adopt a truncated Gaussian prior that is centred on the previous measurement by <cit.> (2010 edition) which is 0.0015. The reported measurement does not include uncertainty, so we use 0.0015  as the standard deviation of the Gaussian prior to allow a broad distribution. We also set a lower bound r_c>1× 10^-4 degree to ensure numerical stability of the calculation of the normalisation factor in the King profile.  For the parameters p_i (and m_i) that specify member fractions in the 4D (1D) model according to Equation <ref>, we adopt a Dirichlet prior over (ϕ_1, …, ϕ_N_ comp) and (c_1,…,c_ N_ comp) which can be interpreted as the uniform distribution under the constraints ∑_i ϕ_i = 1, ϕ_i>0 and c_i = 1, c_i>0.  

Independently for each data set, we obtain random samples from the posteriors using the software package emcee <cit.>, which implements a Markov chain Monte Carlo (MCMC) method.  When fitting the 4D sample, we use 3000 walkers, 5000 steps and thinning parameter of 2.  When fitting the 1D sample, we use 500 walkers, 5000 steps and thinning parameter of 1.
 In both cases, we discard the first half of the chain as burn-in.








§ RESULTS



Random samples sampled from posterior distributions are publicly available at <https://doi.org/10.5281/zenodo.10659516>.  Table <ref> and Table <ref> summarise the marginalised 1D posterior for each model parameter in the 4D and 1D mixture models, respectively, identifying the quantile at (0.1587, 0.5, 0.8413) for each parameter. The member fraction (ϕ_i) for each component in the 4D mixture model, taking the posterior mean, is 95.55% for Sgr, 2.75% for MW, and 1.70% for M54. The member fraction (c_i) for each component in the 1D mixture model, taking the posterior mean, is 71.65% for Sgr, 2.38% for MW, and 25.97% for M54. We note that the tidal radius r_t for M54 tends to be very large in the posterior distribution. We have tried to enforce a strong Gaussian prior using the literature value 0.165 <cit.> as the mean, which will restrict the range of posterior distribution around 0.165, but we do not observe any improvement in the residual distribution or any change in our results. If we take the MAP estimation, the probability mass of normalized King profile within radial distance 8≤ r ≤ 1 is 13%. The distribution is still dominated by the central area. With these considerations, we decide to use the less informative prior which will lead to a large tidal radius. 

We note first that our estimate of the posterior PDF for the Milky Way component of the proper motion distribution is not well converged.  This is due to the relatively small number of MW stars that survive our initial proper motion selection (Figure <ref>; the inferred contribution of MW stars to the 4D model is ≲ 3%).  



Figure <ref> displays the distribution of projected positions in the data (left panel), best-fitting model (middle-left panel), residuals normalised by Poisson uncertainties estimated from the expected values under the best-fitting model (middle-right panel), and the histogram of the normalised residuals (right panel).  Figure <ref> displays the probability distribution for circular radial coordinate R_ circ=√(x^2+y^2), with individual contributions from the modelled M54, Sgr and MW components indicated.  Here the observed data are binned for illustrative purposes only, as the models are fit to the discrete data for stellar position. 
 The lower sub-panel in Figure <ref> shows residuals, normalised by the Poisson uncertainty. 
 We find that our mixture model provides a reasonably good fit, with normalised residuals approximately following a unit-normal distribution (right panel of Figure <ref> and bottom panel of Figure <ref>).  



 





The contours in Figure <ref> represent posterior PDFs that we estimate for spatial centroids (left panel) of M54 and Sgr, and mean proper motions (right panel); for Sgr, for which we allow the mean proper motion to vary linearly with projected positions from the spatial centroid, the mean proper motion is evaluated at the fitted centroid.
The inferred spatial centres of M54 and Sgr are offset by an angle of | Δs|=0.295 ± 0.029^∘—corresponding to a physical separation of 0.135 ± 0.013 kpc at the adopted distance of 26.28 kpc—with M54 projected to the northwest of the Sgr centre.  

The inferred mean (central) proper motions of M54 and Sgr are offset by Δμ= [0.0049, -0.0197] ± [0.0068, 0.0062] mas yr^-1, with M54 moving more slowly than Sgr toward the south.  
However, the precision of mean proper motions estimated for extended objects is limited by spatial covariance within the Gaia's EDR3 proper motion catalog <cit.>.  We use Equation 2 from <cit.> to estimate the magnitude of this systematic error `floor', which depends on the angular size of the object in question.  Using the core radius r_c = 0.0017 for M54 and half-light radius r_h=5.083 for Sgr to represent the relevant angular scales of these objects, the corresponding systematic errors are ϵ_sys≈ 0.026 mas yr^-1 for M54 and ϵ_sys≈ 0.015 mas yr^-1 for Sgr.  Thus the offset that we infer between the mean proper motions of M54 and Sgr is similar to the expected contribution from systematic error. 
Thus the offset that we infer between the mean proper motions of M54 and Sgr has a similar magnitude to the expected systematic error. 








Fig <ref> displays posterior PDFs for mean line-of-sight velocities of M54 and Sgr.  We infer mean values of v_0=139.63 ± 0.92  km s^-1 for M54 and v_0=143.74 ± 0.69  km s^-1 for Sgr (evaluated at the inferred centroid of Sgr), for an offset of Δ v_0=4.1± 1.2 km s^-1.  We infer a significant gradient in Sgr's line-of-sight velocities, e_x=-1.68± 0.57 km s^-1 deg^-1 and e_y=-4.80 ± 0.74 km s^-1 deg^-1, implying a gradient of magnitude ∼ 5 km s^-1 deg^-1 (∼ 10 km s^-1 kpc^-1), similar to the gradient previously reported by <cit.>.

Taken at face value, these results suggest that M54 and Sgr are offset in phase space, with the detected offsets significant at the 10σ and 3.4σ levels in projected sky position and line-of-sight velocity components, respectively.  The inferred proper motion offset, while significant at the ∼ 3σ level based on our posterior PDF alone, is similar to the expected level of systematic error.












§ DISCUSSION





 §.§ Comparison to Previous work


For comparison with our inferences of the centroids and mean motions of M54 and Sgr, Figures <ref> and <ref> indicate previous measurements found in the literature.  Our inference for the spatial centre of M54 is consistent with those reported by <cit.>. However, our inference for the spatial centre of Sgr is offset, by ∼ 15 arcmin, from the center listed in the review by <cit.>, and by similar amounts from measurements reported by <cit.> and <cit.>.  Notably, none of the listed measurements of Sgr's centre agree with each other.  This situation is not entirely surprising, given the different instrumental sensitivities and complex mix of stellar populations near Sgr's centre.  The measurement by <cit.> coincides with a steep `cusp' that they detect at the position of M54, which might reflect contamination of the M giant sample by the relatively young and metal-rich M54 component reported in later work <cit.>.  Indeed <cit.> find that the inclusion of M54 stars in their sample would shift the fitted centroid of Sgr toward M54.  Their measurement that we show in Figure <ref> uses a sample that masks a circular region of radius 0.16^∘ centred on M54, chosen to excise the cluster within its King tidal radius (for comparison, our central mask described in Section <ref> has radius 0.03^∘).  However, <cit.> downplay the resulting offset from M54, noting the potential for systematic errors due to the Gaia scanning pattern that is apparent in their chosen sample from Gaia DR2.  Thus ours is the first inference of an offset between the spatial centres of M54 and Sgr that is based on a mixture model that simultaneously accounts for both objects within the same data set.  Given our selection of a relatively bright magnitude limit (G<17.3; Section <ref>), our data does not show an visible scanning pattern, and thus our result of offset should be less susceptible to systematic errors arising from the spatial dependence of Gaia's faint magnitude limit.

Our inference for M54's proper motion vector is consistent with other recent Gaia-based measurements reported by <cit.>.  However, our inference for Sgr's proper motion (evaluated at Sgr's centre) is statistically inconsistent with those reported by <cit.>, with our measurement indicating a larger component toward the south.  Again we see that all reported measurements for Sgr disagree with each other.  <cit.> note that their PM measurements are affected by Gaia's scanning pattern and varying astrometric incompleteness in Gaia DR2, which contributes to a systematic error that they estimate to be σ_sys∼ 0.030  mas yr^-1 and ∼ 0.036  mas yr^-1 in μ_αcosδ and μ_δ, respectively, similar to the offset between their result and ours.  The measurements by <cit.> and <cit.> are  based on the Gaia DR2 catalogue and the proper motion zero-points are both different and spatially varying between the DR2 and DR3 catalogue and this  systematic error may account for the offset compared to our result.  In any case, even taking any of the available measurements for Sgr's proper motion vector at face value, one finds a similarly significant offset (albeit with reversed direction in the case of the previous measurements based on Gaia DR2) with respect to our inference for the proper motion of M54.  Again, however, we note that the apparent offset is similar in magnitude to the expected contribution from systematic error.

Finally, Fig <ref> shows that previously published measurements of the mean line-of-sight velocities of M54 and Sgr are scattered between the values we infer.  For M54, our estimate (evaluated at the inferred centre of Sgr) is in reasonable agreement with those reported by <cit.>, <cit.> and <cit.>, but more discrepant (with our measurement having a smaller line-of-sight velocity) with those of <cit.> and <cit.>.  For Sgr, our result agrees well with that of <cit.>, marginally with those reported by <cit.> and <cit.>, and poorly (with our measurement being ∼ 4 km s^-1 more positive) with that of <cit.>.  This scatter likely reflects not only systematic errors due to the different zero-points of the different instrumental setups, but also the different criteria used to select and/or identify M54 and Sgr samples.  For example, <cit.> select their Sgr sample from fields observed outside the central region containing M54, with a broad CMD filter chosen to pass a wide range of age and metallicity.  Other studies <cit.> define M54 and/or Sgr samples by applying relatively blue and red filters along the red giant branch sequences; as previously discussed (Section <ref>), this selection can confuse a young, metal-rich component in M54 <cit.> with Sgr, perhaps helping to explain why many of the previous results lie between the mean velocities we measure for the two objects.  Again we emphasise that our mixture model's separation between M54 and Sgr components is based only on the observed distribution of phase-space coordinates, and is agnostic regarding the colour/magnitude and hence age/metallicity properties of each object.  




 §.§ Robustness to sample selection and modelling assumptions


The most significant offset we detect between Sgr and M54 is between their 2D spatial centroids.  We now examine the extent to which this offset may be driven by our choice of modelling assumptions and sample selection.  As described in Section <ref>, we selected our sample of stellar positions (and proper motions) to be spread over a circular field of radius 4^∘ (excluding the crowded region within a field of radius 0.03^∘) centred on M54.  Our model assumes that within this region, the stellar surface density of Sgr is well described by an elliptical Plummer profile.  However, as a result of ongoing tidal interaction with the Milky Way, Sgr is known to display distorted morphology in its outer regions, ultimately extending to tidal tails that encircle the sky <cit.>.  To date, the most detailed view of Sgr's internal stellar structure and kinematics comes from the Gaia-based study by <cit.>, who report a central bar-like feature extending out to the tidal tails, and an overall triaxial morphology, with a slowly rotating (V_ rot∼ 4 km s^-1) central region of radius ∼ 500 pc that is embedded within an outer envelope that appears to be expanding along its longest axis.  

In order to gauge the effect of Sgr's complicated morphology on our inference about its centre, we repeat our fits using samples selected over smaller fields (still centred on M54).  Fig <ref> shows 3σ contours (i.e., enclosing 99.7% of the probability mass) from posterior PDFs for the centroids of M54 and Sgr, obtained using samples spanning fields of radius 4^∘ (the original choice), 3^∘ and 2^∘.  Unsurprisingly, all cases return statistically identical inferences for the centroid of M54.  For Sgr, contours for all cases continue to enclose the centroid inferred from our original sample; however, as the field size decreases, the contours expand and we see a systematic reduction in the inferred offset from M54.  For a field radius of 2^∘, the offset is detected only marginally, with the  3σ contour overlapping the inferred centroid of M54.  


The inflation of statistical uncertainty as the field radius decreases is unsurprising, as the sample becomes more confined to the region well within Sgr's fitted Plummer radius (R_p∼ 5^∘), where the 2D density is close to uniform.  The systematic shift in the inferred Sgr centroid, toward M54 as the field radius decreases, is similar in magnitude to the increase in statistical error.  Thus the systematic component cannot be neglected, and we must acknowledge that the significance with which we detect a spatial offset between M54 and Sgr is vulnerable to model mismatch, particularly with regard to the stellar density distribution toward the outer regions of Sgr.  However, residuals with respect to best-fitting models (Figures <ref> and <ref>) do not reveal obvious signs of model mismatch.  

In order to examine more quantitatively the potential for asymmetric morphology within Sgr to induce spurious detection of an offset from M54, we fit an alternative model for Sgr's 2D spatial configuration that explicitly includes an asymmetric component.

The alternative asymmetric model is constructed by modifying Equation <ref>.  The original matrix R from Equation <ref> is


    R=[       1       0;       0 1/(1-ϵ);         ]R_θ( s-s_0)

Now we construct the asymmetric model to stretch the distribution along just one of the semi-major axes:


    R'=[ |R_00| (1-s_p)/2 s_p + R_00 (1+s_p)/2 s_p;                                      R_10;                                           ]


where R_ij represents the element of the ith row and jth column from R (index starts from 0). This expression ensures that


    R'_00 = R_00     if R_00>=0 
    R_00/s_p      if R_00<0.

We construct the asymmetric model by replacing the R_e^2 = R^⊤R in Equation <ref> with R_e^2 = R'^⊤R'.  The normalisation factor k_P in Equation <ref> is recalculated using the new expression. The new parameter s_p controls how much one side of the Sgr spatial distribution is elongated. The prior on the parameter s_p is U[1, 5].

We fit the asymmetric model using the software package Dynesty <cit.>, which returns an estimate of the marginalised likelihood, or `evidence', as well as a random sample from the posterior.  The evidence provides a metric for model selection.  

The contours in Figure <ref> represent posterior PDFs for the asymmetric stretch factor, s_p, and the inferred centroid of Sgr.  The observed degeneracy indicates that for sufficiently large s_p, near the 2σ contour of our posterior PDF, Sgr's centroid can shift to coincide with the center of M54 at (x,y)=(0,0).  However, comparing results from our asymmetric and original symmetric model, the natural logarithm of the evidence ratio is 3.13± 0.45, in favor of the simpler, symmetric model.  In order to interpret this numerical result, we perform a simulation in which we fit our symmetric and asymmetric models to mock data generated under our best-fitting (to the real data) symmetric model.  In this case we know the input model has no inherent asymmetry, and we obtain a similar evidence ratio, which has a logarithm of 3.03± 0.44, again in favor of the symmetric model.  We conclude that, while the offset we infer between centroids of M54 and Sgr can potentially be caused by asymmetry in Sgr's morphology, the real data do not give reason to prefer at least our simple asymmetric model over the original symmetric one. 












 §.§ Colour/magnitude distribution

We reiterate that our mixture model is agnostic regarding the colour/magnitude distribution of stars belonging to the individual components of M54, Sgr and the Milky Way.  This represents a departure from most previous analyses of M54 and/or Sgr, which typically analyse separate samples that are pre-selected from relatively blue and red sequences defined along the red giant branch <cit.>.  In contrast, our inferences of mixing fraction and structural/kinematic parameters for each mixture component are informed only by the observed stellar positions and motions.  Afterwards, however, we can use the inferred parameters and mixing fractions to examine the corresponding colour/magnitude distribution of each component.  For each of M54, Sgr and the Milky Way, Figure <ref> displays the number of stars per pixel in colour/magnitude space having probability of membership > 50% in that component.  

We see that the Milky Way contribution tends to be fainter and bluer than the M54 and Sgr populations, only sparsely populating the red giant sequences.  Sgr stars lie primarily to the red side of our CMD filter, with a prominent red clump visible at G∼ 16.7.
Finally, we see that while M54 is composed primarily of stars along the blue ridge of our CMD filter, it has a population of high-probability members along the redder sequence as well.  This result is broadly consistent with that of <cit.>, who use MUSE spectroscopy to separate the M54 region into an old (12.2 Gyr), metal-poor ([Fe/H]=-1.41) population a young (2.2 Gyr), metal-rich ([Fe/H]=-0.04) population and an intermediate-age (4.3 Gyr), metal-rich ([Fe/H]=-0.29) population, with the first two being more centrally concentrated than the last.  It is this central concentration that associates the young, metal-rich population with M54 in our model. 
  <cit.> speculate that these stars may have formed in situ in the nuclear region, in response either to the funnelling of gas during Sgr's most recent pericentric passage of the Milky Way, or to the sinking of enriched gas to the centre of a sufficiently massive M54.  





 §.§ Orbits of M54 and Sgr

If M54 is offset from the centre of Sgr in phase space, then it becomes interesting to trace the history of the M54/Sgr interaction, which will depend on the masses and internal structure of both objects as well as the external tidal field.  In order to explore at least some aspects of this dependence, we consider here a simplified `toy' model in which M54 is a point mass and Sgr is embedded within an extended dark matter halo that follows the flexible `CoreNFW' (cNFW) density profile of <cit.>.  This spherically-symmetric halo model is defined by the enclosed-mass profile

    M_ cNFW(r)=f^nM_ NFW(r),

where

    M_ NFW(r)=M_200ln(1+r/r_s)-r/r_s(1+r/r_s)^-1/ln(1+c_200)-c_200(1+c_200)^-1

is the `NFW' profile that characterizes halos produced in cosmological N-body simulations that assume `cold' dark matter <cit.>.  Here, M_200=M_ NFW(r_200) is a proxy for halo mass, with r_200 the radius of the sphere inside which the mean density is 3M_200/(4π r^3_200)=200ρ_c, where ρ_c=3H_0^2/(8π G) is the critical density of the Universe (we assume the Hubble constant has value H_0=70 km s^-1 Mpc^-1).  With r_200, the scale radius r_s defines halo concentration, c_200≡ r_200/r_s.
The function

    f^n≡ [tanh (r/r_c ) ]^n

allows the `cuspy' (dlogρ/dlog r=-1) inner behavior of the NFW profile to be transformed, to a degree controlled by parameter 0<n≤ 1, to a `core' of uniform density within a sphere of radius r_c.

In order to examine sensitivity to the structure of Sgr's dark matter halo, we consider two possibilities, both of which assume that M54 is a point particle of mass M_ M54=1.41× 10^6 M_⊙ <cit.>  and Sgr consists of a stellar component, described by a spherical Plummer model, embedded within a cNFW dark matter halo.
In the first case, the cNFW halo has r_c=0, which simplifies to the original NFW profile.  In the second case, the cNFW halo has a fully-formed (n=1) core of radius, r_c=r_p, that is assumed to equal the Plummer radius fit to the observed stellar density profile <cit.>.
In both cases, we adopt a dark matter halo concentration of c_200=10, and for Sgr a total stellar mass (within the bound Plummer sphere) of M_*, Sgr=2× 10^7 M_⊙ <cit.>.
For both dark matter profiles, we choose M_200 so that the dark matter mass enclosed within the Plummer radius matches dynamical mass estimated as M(r_half) = 2.5 r_pσ_v^2/G <cit.>, where we adopt r_p=5.092^∘ as the MAP result from our modeling and σ_v=11.4 km s^-1 is the velocity dispersion <cit.>.  For the NFW and cNFW profiles, this requirement gives M_200=1.29 × 10^9  M_⊙ and M_200=2.10 × 10^9  M_⊙, respectively.  
In order to account for any effect of the Large Magellanic Cloud (LMC) on the orbits of M54 and Sgr, we include an LMC potential that we model using a Hernquist profile with total mass 1.38 × 10^11 M_⊙ and scale radius 16.09  kpc
<cit.>; however, we find that the inclusion or exclusion of the LMC potential has negligible effect on our results.  

We sample the present-day 2D position, proper motions and line-of-sight velocities for both M54 and Sgr from the posterior PDFs inferred from our mixture models.
We include the systematic errors <cit.> in the proper motion samples by adding Gaussian noise to the proper motion components sampled from our posterior PDF, accounting for uncertainty `floors' of ϵ_sys= 0.026 mas yr^-1 for the proper motion components of M54 and ϵ_sys=0.015 mas yr^-1 for Sgr.
We draw present-day line-of-sight distances to M54 and Sgr independently from a normal distribution with mean 26.28 kpc and standard deviation 0.33 kpc, corresponding to the measurement by <cit.> of the distance to M54.  
For Sgr, assumed here to be spherically symmetric, we sample the Plummer radius from our posterior PDF for r_h.  For the LMC, we adopt the present-day phase space coordinates from <cit.>.


Given the stated initial conditions, we use the software package  <cit.> to calculate the orbits of M54, Sgr and the LMC within a Milky Way potential that adopts 's MilkyWayPotential class, which assumes a disk model from <cit.>.  We choose 's RK5Integrator, which uses a fifth order Runge-Kutta method.  In order to gauge the effects of dynamical friction on M54 from both Sgr and the MW and on Sgr and the LMC from the MW, we calculate orbits with and without using galpy's ChandrasekharDynamicalFrictionForce class <cit.>, which implements the Chandrasekhar dynamical friction force following <cit.>.  For each combination of NFW vs. cNFW for Sgr's dark matter profile and with/without dynamical friction implemented, we draw 100 samples from our posterior PDFs that provide initial conditions, and then calculate the corresponding orbits over the past 1 Gyr. The timestep is 5 Myr and the number of steps is 200.  

For each computed orbit, Fig <ref> depicts the recent history of the (3D) spatial offset between M54 and the centre of Sgr, tracing back 1 Gyr from the present day.  Left and right columns show results using the NFW and coreNFW model, respectively, for Sgr's dark matter halo; top and bottom panels show results obtained without and with dynamical friction implemented.

For cases where the present-day offset is ≲ 0.7 kpc, we observe a dependence on the dark matter density profile within Sgr.  If Sgr's dark matter halo follows the NFW profile, we find that most orbits are approaching apocentre, with the M54/Sgr separation increasing from an offset as small as ∼ 0.1 kpc at the most recent pericentric passage, which occurred between ∼ 15 - 40 Myr ago.  In the absence of dynamical friction, the orbit remains approximately unchanged for the past Gyr for most of the samples (upper left panel of Figure <ref>).  With dynamical friction implemented (lower left panel of <ref>), the offset at apocentre has decreased steadily from a value of ∼ 1.5 kpc over the past 700 Myr. 

In contrast, if Sgr's dark matter halo follows the adopted coreNFW profile, then present-day offsets that are ≲ 0.7 kpc are close to pericentre for most of the samples, after a most recent apocentre ≳ 100 Myr ago, at which the offset was ≳ 0.5 kpc (right panels of Figure <ref>).  As expected, the implementation of dynamical friction has relatively little effect when the dark matter halo is cored <cit.>, allowing M54 to `wobble' about the centre indefinitely <cit.>.

If the present-day offset is ≳ 0.7 kpc, and thus dominated by a large offset in line-of-sight distance between M54 and Sgr, all sampled orbits would imply that M54 fell into Sgr within the past 200 Myr, regardless of whether Sgr's dark matter halo follows a NFW or cNFW profile.  However, this scenario seems unlikely, as one expects Sgr to be losing mass, not accreting, during this timeframe that includes Sgr's own pericentric passage about the Milky Way <cit.>.

We emphasise that our orbit calculations neglect realistic details of the Sgr/M54 system—most notably, Sgr's tidal disruption by the Milky Way.  The orbits that we calculate imply that Sgr's most recent pericentric passage about the Milky Way brought it within ∼ 15 kpc of the Galactic centre, just ∼ 45 Myr ago, consistent with the Sgr orbit inferred by <cit.>.  However, even if we restrict our calculations only to the previous 50 Myr, the qualitative conclusion remains that an M54/Sgr system that has internal offset ≲ 0.7 kpc tends to be appropaching apocentre if Sgr is embedded within an NFW halo, and near pericentre if Sgr is embedded within a cored halo.  A more thorough investigation of the history of the M54/Sgr interaction would require updating previous N-body approaches <cit.> in light of the present-day phase-space coordinates presented above.








§ SUMMARY & CONCLUSIONS



Motivated by longstanding questions regarding the formation and orbital histories of nuclear star clusters, here we have analyzed the 5D (sky positions and proper motions from Gaia EDR3, line-of-sight velocities from APOGEE DR17) structure and kinematics of the M54/Sgr system, the only example of a NSC/host system for which such multi-dimensional data are available.  In order to obtain simultaneous estimates of the centers of M54 and Sgr in each of these coordinates, we have constructed and fit a mixture model that accounts also for contamination from the Milky Way foreground/background.  Our approach differs from most previous efforts, which consider M54 and/or Sgr samples separately after filtering according to colour/magnitude and/or 2D spatial criteria.  In analysing both M54 and Sgr populations simultaneously within the same sample, our inferences about offsets between the two objects incur fewer systematic errors that arise when comparing estimates made using either independent data sets or independent subsamples of those data sets.  


We have reported statistically significant offsets in the spatial centres and line-of-sight velocities of M54 and Sgr. The offset between the spatial centres is 0.295± 0.029 degrees, corresponding to a projected (2D) offset of 0.135± 0.013 kpc at the distance of 26.28 kpc.  The offset between the line-of-sight velocities is 4.1± 1.2 km s^-1.  The offset between the PM centres of the M54 and Sgr dSph is [Δμ_αcosδ,Δμ_δ]=[0.0049, -0.0197] ± [0.0068, 0.0062] mas yr^-1, comparable to the expected systematic error.  

By far the most significant offset we infer is between the spatial centres of M54 and Sgr.  However, this result comes with the caveat that the degree of statistical significance is sensitive to our choice of sample field size.  This sensitivity points to the possibility that the inferred offset may be driven by large-scale asymmetry of the spatial distribution of stars within Sgr.  While our fiducial model allows flattened elliptical morphology, we find that the inferred offset between M54 and Sgr can be reduced if we further allow for bilateral asymmetry.  However, we find that standard criteria for model selection disfavour our asymmetric model over the fiducial one.  We conclude, then, that the present data are reasonably well fit by our fiducial model, with M54 significantly offset from Sgr, and do not favour an alternative in which M54 is located exactly at the centre of an asymmetric Sgr host.

Based on our estimates of the 5D phase-space coordinates of M54 and Sgr at the present day, we also investigated the recent orbital history of the interacting pair.  Treating M54 as a point mass and neglecting mass lost from Sgr due to tidal disruption, we calculate orbits within a potential that includes both objects as well as the Milky Way and LMC.  We find that if Sgr is embedded in a cuspy (NFW) dark matter halo, then the M54/Sgr orbit is currently most likely approaching apocentric when the present-day offset is ≲ 0.7 kpc, having decayed under dynamical friction from earlier apocentres ≳ 1 kpc.
In contrast, if Sgr's dark matter halo has a central core, then the M54/Sgr orbit is near pericentre if the present-day offset is ≲ 0.7 kpc.
If the present-day offset between M54 and Sgr is ≳ 0.7 kpc, and thus dominated by an offset in the line-of-sight distance, the calculated orbits would imply that M54 fell into Sgr within the past 200 Myr regardless of whether Sgr's dark matter halo is cored or cusped.

Finally, the results of our mixture model imply that in addition to a dominant population of old, metal-poor stars, the M54 component includes a population of young, metal-rich stars.  This finding is consistent with the MUSE spectroscopy of <cit.>, who associate a young, metal-rich and spatially compact stellar population with M54, and supports the suggestion by <cit.> that these stars may have formed in response to the sinking of enriched gas into the center of a massive M54, perhaps triggered during Sgr's most recent pericentric passage of the Milky Way.  Further exploration of this scenario in future work will require N-body simulations, which can be constrained by the results presented above.





§ ACKNOWLEDGEMENTS

This work is supported by National Science Foundation (NSF) grants AST-1909584 and AST-2206046.  This work made use of the Q3C software <cit.>, NumPy <cit.>, SciPy <cit.>, emcee <cit.>, Numba <cit.>, gala <cit.>, Dynesty <cit.>, galpy [http://github.com/jobovy/galpy] <cit.> and Astropy[http://www.astropy.org] which is a community-developed core Python package and an ecosystem of tools and resources for astronomy <cit.>.

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.




§ DATA AVAILABILITY

The data underlying this article are available in <https://doi.org/10.5281/zenodo.10659516>. The datasets were derived from the Gaia EDR3 and the APOGEE DR17. The the Gaia EDR3 is publicly available at <https://www.cosmos.esa.int/web/gaia/earlydr3>, and the APOGEE DR17 is publicly available at <https://www.sdss4.org/dr17/>.



mnras
 










