Draft version April 23, 2021

Typeset using LATEX twocolumn style in AASTeX63

The Age–Metallicity–Specific Orbital Energy Relation for the Milky Way’s Globular Cluster System
Confirms the Importance of Accretion for Its Formation
Turner Woody1 and Kevin C. Schlaufman1
1 Departmentof Physics and Astronomy
Johns Hopkins University
arXiv:2104.10697v1 [astro-ph.GA] 21 Apr 2021

3400 N Charles St
Baltimore, MD 21218, USA

(Received October 23, 2020; Revised March 23, 2021; Accepted April 11, 2021)

Accepted for publication in the Astronomical Journal

ABSTRACT
Globular clusters can form inside their host galaxies at high redshift when gas densities were higher
and gas-rich mergers were common. They can also form inside lower-mass galaxies that have since
been accreted and tidally disrupted, leaving their globular cluster complement bound to higher-mass
halos. We argue that the age–metallicity–specific orbital energy relation in a galaxy’s globular cluster
system can be used to identify its origin. Gas-rich mergers should produce tightly bound systems in
which metal-rich clusters are younger than metal-poor clusters. Globular clusters formed in massive
disks and then scattered into a halo should have no relationship between age and specific orbital
energy. Accreted globular clusters should produce weakly bound systems in which age and metallicity
are correlated with each other but inversely correlated with specific orbital energy. We use precise
relative ages, self-consistent metallicities, and space-based proper motion-informed orbits to show that
the Milky Way’s metal-poor globular cluster system lies in a plane in age–metallicity–specific orbital
energy space. We find that relatively young or metal-poor globular clusters are weakly bound to
the Milky Way, while relatively old or metal-rich globular clusters are tightly bound to the Galaxy.
While metal-rich globular clusters may be formed either in situ or ex situ, our results suggest that
metal-poor clusters formed outside of the Milky Way in now-disrupted dwarf galaxies. We predict that
this relationship between age, metallicity, and specific orbital energy in a L∗ galaxy’s globular cluster
system is a natural outcome of galaxy formation in a ΛCDM universe.

Keywords: Galaxy formation(595) — Globular star clusters(656) — Milky Way dynamics(1051) —

Milky Way Galaxy(1054) — Milky Way formation(1053) — Milky Way stellar halo(1060)

1. INTRODUCTION ferent from its current z = 0 state. Indeed, the in-

Globular cluster formation is one the enduring un- terval spanning the 16th and 84th percentiles of the
solved problems in astrophysics. Globular clusters can Milky Way globular cluster system’s age distribution is
form either in situ inside their parent galaxy or ex situ in 11 Gyr . τ . 13 Gyr or 2.4 . z . 7.3 in terms of red-
lower-mass galaxies that have been accreted and tidally shift (e.g., Marı́n-Franch et al. 2009; VandenBerg et al.
disrupted by a more massive halo (e.g., West et al. 2013; Wright 2006). Both in situ and ex situ formation
2004; Forbes et al. 2018). In either case, the ancient channels are theoretically expected to contribute to the
ages inferred for most globular clusters indicate that Milky Way’s globular cluster system (e.g., Griffen et al.
they formed long ago when the universe was very dif- 2010; Renaud et al. 2017).
In the in situ formation scenario, globular clusters
form inside their parent galaxy in gas-rich galaxy merg-
Corresponding author: Turner Woody ers (e.g., Ashman & Zepf 1992; Muratov & Gnedin 2010;
[email protected] Li & Gnedin 2014, 2019; Kim et al. 2018) and/or in
gas-rich disks (e.g., Kravtsov & Gnedin 2005; Kruijssen
2 Woody & Schlaufman

2015; Pfeffer et al. 2018; Keller et al. 2020). Galax- rich and metal-poor globular clusters can form in mas-
ies in general and galaxies destined to be similar to the sive dSph galaxies, low-mass dSph galaxies will only be
Milky Way in particular are observed to be smaller at able to form metal-poor clusters. More sophisticated ar-
high redshift when globular clusters form in this scenario guments based on cosmological dark matter-only or hy-
(e.g., Patel et al. 2013a,b). Consequently, globular clus- drodynamical simulations have reached similar conclu-
ters formed in this way should be tightly bound to their sions (e.g., Choksi et al. 2018; Kruijssen 2019; Kruijssen
now more-massive and larger parent galaxies at z = 0 et al. 2019a,b).
(e.g., Leaman et al. 2013). Another robust prediction The mass of an accreted globular cluster’s parent
of the gas-rich merger scenario is that metal-rich clus- galaxy will also affect its orbit. A dwarf galaxy orbiting
ters should be systematically younger than metal-poor inside the dark matter halo of a more massive galaxy will
clusters (Li & Gnedin 2014, 2019). lose angular momentum due to dynamical friction and
In the ex situ scenario, globular clusters form outside its orbit will decay on a timescale that is proportional
their z = 0 host halo in parent galaxies that have since to its mass (e.g., Binney & Tremaine 2008). Massive
been accreted and tidally disrupted (e.g., Searle & Zinn dSph galaxies will therefore be dragged into the inner
1978; Mackey & Gilmore 2004; Forbes & Bridges 2010; regions of their host halo and tidally disrupted much
Forbes 2020; Massari et al. 2019). While in most cases more quickly than low-mass dSph galaxies. In partic-
it is difficult to conclusively associate a candidate ac- ular, a massive dSph galaxy resembling Fornax with a
creted globular cluster with its parent galaxy, the glob- total mass Mtot ∼ 108 M capable of forming a globu-
ular clusters Arp 2, NGC 6715 (M 54), Pal 12, Terzan 7, lar cluster with [Fe/H] ≈ −1.0 will be tidally disrupted
and Terzan 8 have all been securely associated with the 100 times faster than a low-mass dSph galaxy with a
Sagittarius dwarf spheroidal (dSph) galaxy (e.g., Law & total mass Mtot ∼ 106 M that could not form a glob-
Majewski 2010; Sohn et al. 2018). Given the collision- ular cluster more metal-rich than [Fe/H] ≈ −2.0 (e.g.,
less nature of their accretion, globular clusters formed Pascale et al. 2018; Simon 2019). The metal-rich globu-
ex situ should be on average less tightly bound to their lar clusters deposited deep inside the more massive halo
host galaxy at z = 0 than clusters formed in situ in when their parent massive dSph galaxy is tidally dis-
gas-rich conditions. rupted will be much more tightly bound to the massive
The observed properties of accreted globular clusters halo than any metal-poor globular clusters left behind
can be used to set limits on the properties of their parent when their parent low-mass dSph galaxy is tidally dis-
galaxies. Likewise, the ensemble properties of a galaxy’s rupted much later.
accreted globular cluster system can be used to explore The ages of accreted globular clusters should be re-
its formation. The maximum metallicity realized in a lated to the masses and therefore the metallicities of
dwarf galaxy’s globular cluster system should not ex- their parent galaxies too. In similar environments,
ceed the dwarf galaxy’s typical metallicity. This is the galaxies destined to be massive form in higher σ peaks in
case in both the Fornax and Sagittarius dSph galax- the primordial matter density distribution than galaxies
ies. The Fornax dSph has h[Fe/H]i = −1.04 (Kirby destined to be low mass (e.g., Mo et al. 2010). Since the
et al. 2013), while its five globular clusters span the dynamical time tdyn scales like tdyn ∝ ρ−1/2 , in similar
range −2.5 . [Fe/H] . −1.4 (e.g., Letarte et al. 2006; environments more massive galaxies form stars earlier
Larsen et al. 2012). While the Sagittarius dSph has and therefore the ages of the oldest stellar populations
h[Fe/H]i ≈ −0.5 (e.g., Chou et al. 2007; Hasselquist et al. in galaxies should be correlated with their masses. The
2017), the five globular clusters securely associated with implication is that in similar environments the age of
it span the metallicity range −2.3 . [Fe/H] . −0.6 (e.g., the oldest globular cluster formed in a dwarf galaxy
Cohen 2004; Sbordone et al. 2007; Mottini et al. 2008; should be correlated with its mass. This is the case
Carretta et al. 2010, 2014). A lower limit on the stellar for the Mtot ≈ 9 × 107 M Fornax dSph galaxy, the
mass of a galaxy necessary for the formation of a glob- Mtot & 4 × 108 M Sagittarius dSph galaxy, and the
ular cluster of a given metallicity can therefore be de- Mtot ∼ 1011 M Large Magellanic Cloud (e.g., Pas-
rived from a stellar mass–metallicity relationship. Kirby cale et al. 2018; Vasiliev & Belokurov 2020; Erkal et al.
et al. (2013) showed that there exists a linear relation- 2019). The Fornax globular clusters 1, 3, and 5 are com-
ship between log stellar mass and metallicity for dwarf parable in age to NGC 4590 (M 68) while cluster 4 is
galaxies over five orders of magnitude such that the typ- somewhat younger (e.g., Buonanno et al. 1998, 1999).
ical metallicity of a dwarf galaxy rises linearly with log On the other hand, Terzan 8 is the oldest globular clus-
stellar mass from [Fe/H] ≈ −2.5 at M∗ ∼ 103 M to ter with a confirmed Sagittarius association and is older
[Fe/H] ≈ −1.0 at M∗ ∼ 108 M . So while both metal- than NGC 4590 (e.g., Marı́n-Franch et al. 2009; Sohn
 The Age–Metallicity–Specific Orbital Energy Relation for MW Globular Clusters 3

et al. 2018). Likewise, the oldest of the Large Magellanic 2018; Lindegren et al. 2018; Luri et al. 2018). We sup-
Cloud’s globular clusters that are not affected by crowd- plement the globular cluster proper motions published
ing, Hodge 11, is older than Terzan 8 but younger than in Gaia Collaboration et al. (2018a) with the Gaia DR2-
the Milky Way’s oldest globular clusters (Wagner-Kaiser based proper motions for NGC 6584 and NGC 6723 from
et al. 2017). Analyses based on scaling relations cal- Baumgardt et al. (2019) that have proper motion pre-
ibrated with cosmological hydrodynamical simulations cisions µ/σµ > 85.1 We also add proper motions mea-
reach similar conclusions (e.g., Kruijssen 2019). surements for NGC 1261, NGC 4147, NGC 6101, Terzan
The relationship between age and metallicity in the 7, Arp 2, Terzan 8, NGC 6934, and Pal 12 from the
Galaxy’s globular cluster system has long been used to High-resolution Space Telescope Proper Motion Collab-
constrain Milky Way formation. As we described above, oration (HSTPROMO; Sohn et al. 2018). We list our
analyses of the relationship between age, metallicity, and input globular cluster proper motions in Table 1.
orbital properties in the Milky Way’s globular cluster Globular cluster age determinations are inherently dif-
system can be even more informative. If globular clus- ficult to put on an absolute age scale and can be non-
ters formed in situ are produced in gas-rich mergers, trivially affected by even small differences in analysis
then metal-rich clusters should be younger than metal- techniques. Since the analyses we will describe in Sec-
poor clusters and both metal-rich and metal-poor clus- tion 3 only make use of relative age differences, we
ters should be relatively tightly bound to their parent choose to use the homogeneous relative globular cluster
galaxy. If globular clusters are formed in situ in gas- ages presented in Marı́n-Franch et al. (2009) to avoid
rich disks and subsequently scattered into the halo, then the difficulties associated with absolute ages. Marı́n-
there should be no relationship between age and specific Franch et al. (2009) calculated relative clusters ages
orbital energy. In the ex situ formation and accretion τ = age/12.8 Gyr for two different metallicity scales
scenario, older, metal-rich globular clusters should be (Zinn & West 1984; Carretta & Gratton 1997) as re-
tightly bound to their z = 0 host galaxies while younger, ported in Rutledge et al. (1997a,b) and four different
metal-poor globular clusters should be weakly bound to sets of theoretical stellar models (Bertelli et al. 1994;
their z = 0 host. These straightforward arguments have Girardi et al. 2000; Pietrinferni et al. 2004; Dotter et al.
been confirmed by cosmological hydrodynamical simu- 2007). We list our adopted globular cluster metallici-
lations (e.g., Pfeffer et al. 2020). ties and normalized ages in Table 2. While we focus on
In this paper, we calculate the orbits of the Milky the relative ages calculated assuming the more modern
Way’s globular clusters with precise space-based astrom- Carretta & Gratton (1997) metallicity scale and Dot-
etry. We use those orbits along with precise relative ages ter et al. (2007) models given in Table 2, the analysis
and metallicities from the literature to quantify the re- we describe in Section 3 reaches the same conclusion re-
lationship between age, metallicity, and specific orbital gardless of the assumed metallicity scale or model grid.
energy in the Milky Way’s globular cluster system. We We adopt the globular cluster distances and radial ve-
describe the construction of our globular cluster sample locities listed in Table 1 from the December 2010 revi-
in Section 2. We detail in Section 3 our orbital inte- sion of the Harris (1996) compilation.2 Since the De-
grations and statistical analyses of the age–metallicity– cember 2010 revision of the Harris (1996) compilation
specific orbital energy relation in our complete sample does not provide uncertainties for its distance estimates,
of globular clusters as well as its metal-rich/metal-poor we estimate distance uncertainties by first transforming
subsamples. We review our results and their implica- the distance inferences into distance moduli. We next
tions for the formation of the Milky Way’s globular clus- assume a uniform distance modulus uncertainty of 0.05
ter system specifically and globular cluster formation in mag and sample each cluster’s distance modulus from
general in Section 4. We conclude by summarizing our a normal distribution centered on its distance modu-
findings in Section 5. lus with standard deviation 0.05 mag. We then trans-
form the simulated distance moduli calculated in this
2. DATA way back into distances and uncertainties. We will use
the 57 Milky Way globular clusters listed in Tables 1 and
In an analysis of the age–metallicity–specific orbital
2 with homogeneous relative ages and metallicities plus
energy relation in the Milky Way’s globular cluster sys-
tem, the limiting factors are precise cluster ages and
1
proper motions. Most of our proper motion measure- We use this threshold because all of the globular cluster proper
ments come from Gaia Data Release 2 (DR2; Gaia motions from Gaia Collaboration et al. (2018a) have µ/σµ > 85
in the Baumgardt et al. (2019) catalog.
Collaboration et al. 2016, 2018b,a; Crowley et al. 2016; 2 https://physwww.mcmaster.ca/∼harris/mwgc.ref
Fabricius et al. 2016; Arenou et al. 2018; Hambly et al.
4 Woody & Schlaufman

precise space-based proper motions to explore the re- disk is represented by a Miyamoto–Nagai potential with
lationship between age, metallicity, and specific orbital a radial scale length of 3 kpc and a vertical scale height of
energy in Section 3. 280 pc (Miyamoto & Nagai 1975). Its halo is modeled as
a Navarro–Frenk–White (NFW) halo with a scale length
3. ANALYSIS
of 16 kpc (Navarro et al. 1996). The scaled version of
the MWPotential2014 potential differs from the default
We use the right ascensions α, declinations δ, proper MWPotential2014 potential in that we have increased
motions µα & µδ , distances d, and radial velocities vr its halo mass by 50%. This change increases the virial
given for the 57 globular clusters listed in Table 1 to mass to Mvir = 1.2 × 1012 M . In both cases, we set the
calculate their orbits and specific orbital energies. We solar distance to the Galactic center to R0 = 8.122 kpc,
use a Monte Carlo approach to average over the un- the circular velocity at the Sun to V0 = 238 km s−1 ,
certainties in our input proper motions, distances, and and the height of the Sun above the plane to z0 = 25
radial velocities. On each iteration, we sample the nec- pc (Jurić et al. 2008; Bland-Hawthorn & Gerhard 2016;
essary input data for an orbit integration (α, δ, µα , µδ , Gravity Collaboration et al. 2018).
d, and vr ) from the data and uncertainties in Table 1. The McMillan17 potential has Mvir = 1.4 × 1012 M .
When possible, we account for the covariance between Its bulge is parameterized as a spherical power law den-
simultaneously inferred proper motion and parallax in- sity profile with exponent −1.8 that is finite at r = 0
ferences. For the globular clusters with proper motions and has an exponential cut-off at 2.1 kpc. Its disk is
from Gaia Collaboration et al. (2018a), we account for split into four components: a thin stellar disk, a thick
the covariance between proper motions and parallaxes stellar disk, an H I gas disk, and an H2 gas disk. Its
even though we use Harris (1996)-based distances in- two stellar disks are modeled as exponential in both the
stead of the Gaia DR2 parallax-based distances. For the cylindrical radius R and height z, with scale lengths
globular clusters with proper motions from Baumgardt Rd,thin = 2.6 kpc, Rd,thick = 3.6 kpc, zd,thin = 300 pc,
et al. (2019), we account for the covariance between and zd,thick = 900 pc. Its two gas disks are modeled
the proper motion components. Sohn et al. (2018) did as exponential in R with a central hole and a modified
not provide estimates of the covariance between proper exponential profile in the z coordinate. Its halo is mod-
motion components, so we sample each proper motion eled as a NFW halo with a scale length of 19.6 kpc. The
component independently from within its quoted uncer- McMillan17 potential uses a solar distance and circular
tainty. For distances and radial velocities, we sample velocity of 8.21 kpc and 233.1 km s−1 .
from normal distributions with the means and standard The Milky Way’s mass within 25 kpc M25 is especially
deviations given in Table 1. We repeat this process 100 relevant for this problem, and estimates of M25 for the
times to generate 100 sets of initial conditions for each Milky Way based on globular cluster kinematics suggest
globular cluster. that M25 ≈ 0.26+0.10 12
−0.06 × 10 M (Eadie & Jurić 2019).
For each Monte Carlo initial condition realization, we The M25 values for the default MWPotential2014, the
integrate a cluster’s orbit in model Milky Way poten- scaled MWPotential2014, and the McMillan17 poten-
tials using the galpy python module3 (Bovy 2015). We tials are M25 = 0.25 × 1012 M , M25 = 0.29 × 1012 M ,
calculate orbits from the present five Gyr into the past and M25 = 0.28 × 1012 M respectively. All three
with a time step of 200,000 years assuming a solar mo- potentials are therefore consistent the best available
tion with respect to the local standard of rest (U , V , Milky Way constraints in the same volume. While both
W ) = (11.1, 12.24, 7.25) km s−1 (e.g., Schönrich et al. the default MWPotential2014 and McMillan17 poten-
2010). tials are generally consistent with literature estimates
We calculate orbits in three different model Milky Way of the Milky Way’s virial mass, most recent mass esti-
potentials: the default MWPotential2014 from Bovy mates making use of Gaia DR2 proper motions suggest
(2015), a scaled version of MWPotential2014 with its virial masses in the range 1.1 × 1012 M . Mvir .
halo mass increased by 50%, and McMillan17 from 1.4 × 1012 M (e.g., Wang et al. 2020). While the virial
McMillan (2017). In the default MWPotential2014 po- mass of the McMillan17 potential is in this interval,
tential, the virial mass and radius are Mvir = 0.8 × the virial mass of the default MWPotential2014 poten-
1012 M and rvir = 245 kpc respectively. Its bulge is tial is below this range. In addition, the McMillan17
parameterized as a power-law density profile with expo- potential provides a more comprehensive model of the
nent −1.8 that is exponentially cut-off at 1.9 kpc. Its Milky Way’s disk than either the default or scaled
MWPotential2014 potentials. This difference could be
3 https://github.com/jobovy/galpy important for the relatively nearby clusters that dom-
 The Age–Metallicity–Specific Orbital Energy Relation for MW Globular Clusters 5

inate our input sample. For all of these reasons, the We use F tests to compare the full age–metallicity–
McMillan17 potential likely provides a better approxi- specific orbital energy model we derived above with
mation to the true Milky Way potential than the default each of its nested submodels to verify that the full age–
and scaled MWPotential2014 potentials. metallicity–specific orbital energy model best describes
Our orbit integrations produce 100 plausible orbits for the data in Table 2. For each subsample in each poten-
each of the 57 globular clusters listed in Table 1 in the tial, we tested three null hypotheses: βτ = 0, βM = 0,
default MWPotential2014, scaled MWPotential2014, and βτ = βM = 0. We report the results of those model
and McMillan17 potentials. We calculate the 16th, comparisons in Table 4. For the metal-poor subsam-
50th, and 84th percentiles of the specific orbital energy ple in all three potentials, the p values for each test are
(SOE) distribution for every cluster and use those data small enough to confidently reject the null hypotheses.
to define each cluster’s specific orbital energy and un- The implication is that the full linear model is a better
certainty in each potential. We list the specific orbital description of the metal-poor subsample than any of its
energies derived in this way for each cluster in Table 2. nested submodels. In the case of the McMillan17 po-
We use the statsmodels python module (Seabold & tential, the full linear model is a better description of
Perktold 2010) to fit an ordinary least squares linear re- the complete sample than any of its nested submodels.
gression of the form SOE = β0 + βτ τ + βM [M/H] to On the other hand, for the metal-rich subsample there
the age–metallicity–specific orbital energy distribution is no reason to prefer the full age–metallicity–specific
produced by each set of initial conditions in all three orbital energy model over any of its nested submodels.
potentials. We aggregate the results and calculate the We also use the Bayesian information criterion (BIC)
16th, 50th, and 84 percentiles of each distribution. We to compare models. Smaller BIC values indicate better
report the results in Table 3. In all three potentials, models, and the results in Table 5 show that the full age–
we find evidence of an age–metallicity–specific orbital metallicity–specific orbital energy model is preferred to
energy relationship in which relatively old, metal-rich any of its nested submodels for both the complete sam-
globular clusters have more negative specific orbital en- ple and metal-poor subsample in the McMillan17 poten-
ergies and are therefore more tightly bound to the Milky tial. All models are approximately similar in the default
Way than relatively metal-poor, young globular clusters. and scaled MWPotential2014 potentials.
This relationship is especially strong in the McMillan17 We plot the distribution of the Milky Way’s globu-
potential, as βτ and βM in that potential are significant lar cluster systems in the age–metallicity–specific or-
at approximately the 5-σ level. bital energy space for all three potentials in Figure 1.
We separately study the age–metallicity–specific or- It is clear that the 57 globular clusters listed in Ta-
bital energy relationship for metal-rich and metal-poor bles 1 and 2 lie in a plane in age–metallicity–specific
clusters. We split the complete sample into metal- orbital energy space in the McMillan17 potential, ex-
rich/metal-poor subsamples at [M/H] = −0.8 because actly as predicted by our linear modeling. We plot
that cleanly separates the two peaks of the complete two-dimensional projections of age–metallicity–specific
sample’s bimodal metallicity distribution. We find sim- orbital energy space for the default MWPotential2014,
ilar results for both the complete sample and the metal- scaled MWPotential2014, and McMillan17 potentials in
poor subsample of globular clusters that we define as Figures 2, 3, and 4. There remains significant disper-
the 45 globular clusters with [M/H] ≤ −0.8. The t sion about the linear model visualized in Figures 2, 3,
values and therefore the statistical significances of the and 4, so the age–metallicity–specific orbital energy re-
coefficients defining the relation are generally larger in lation inference we described in Section 1 does not fully
absolute value for the metal-poor subsample. For the explain the highly nonlinear process of globular cluster
12 metal-rich clusters with [M/H] > −0.8, the age– formation. Nevertheless, the meaningfully non-zero R2
metallicity–specific orbital energy relationship is much values presented in Table 3 show that a simple three-
weaker. Indeed, if we exclude the two metal-rich glob- parameter linear model explains at least 15%—and pos-
ular clusters associated with the Sagittarius dSph, then sibly as much as 38%—of the variance in the metal-poor
in both the default and scaled MWPotential2014 poten- globular cluster age–metallicity–specific orbital energy
tials the coefficients defining the age–metallicity–specific relation.
orbital relation energy change signs such that older,
more metal-rich clusters are less tightly bound to the 4. DISCUSSION
Milky Way. In this case, there is no significant age–
We find a significant age–metallicity–specific orbital
metallicity–specific orbital energy relationship for the
energy relationship for the metal-poor subsample of 45
the McMillan17 potential.
globular clusters with [M/H] ≤ −0.8 listed in Tables
6 Woody & Schlaufman

Table 3. Median Regression Model Coefficients and t Values

Sample β0 β0 βτ βτ βM βM R2

(104 ) t Value (104 ) t Value (104 ) t Value

Default MWPotential2014
Complete -5.39 ± 2.32 -2.33 -4.12 ± 2.54 -1.62 -1.16 ± 0.63 -1.84 0.081
Metal-poor -2.21 ± 3.07 -0.72 -8.64 ± 3.62 -2.39 -1.98 ± 0.92 -2.16 0.146
Metal-richa -11.26 ± 5.19 -2.17 -1.38 ± 4.20 -0.33 -7.45 ± 8.09 -0.92 0.087
Metal-richb -22.52 ± 10.02 -2.25 19.30 ± 10.80 1.79 9.34 ± 10.16 0.92 0.317
Scaled MWPotential2014
Complete -10.17 ± 2.44 -4.16 -3.95 ± 2.68 -1.47 -1.40 ± 0.67 -2.10 0.089
Metal-poor -5.89 ± 3.19 -1.85 -9.53 ± 3.76 -2.54 -2.15 ± 0.95 -2.26 0.160
Metal-richa -18.20 ± 5.42 -3.36 -0.12 ± 4.38 -0.03 -9.53 ± 8.46 -1.13 0.131
Metal-richb -30.52 ± 9.72 -3.14 22.91 ± 10.48 2.19 9.52 ± 9.85 0.97 0.407
McMillan17
Complete 1.39 ± 3.95 0.35 -24.52 ± 4.33 -5.66 -4.85 ± 1.08 -4.49 0.435
Metal-poor 2.80 ± 5.23 0.54 -28.18 ± 6.17 -4.57 -6.32 ± 1.56 -4.05 0.381
Metal-richa -3.86 ± 9.24 -0.42 -24.76 ± 7.47 -3.32 -15.79 ± 14.41 -1.10 0.551
Metal-richb -38.61 ± 19.16 -2.02 15.39 ± 20.67 0.745 -5.92 ± 19.43 -0.31 0.123
Note—Median age–metallicity–specific orbital energy model coefficients, coefficient uncertainties, t values, and R2
values calculated from 100 Monte Carlo iterations for the complete, metal-poor, and metal-rich subsamples (both
including and excluding clusters associated with the Sagittarius dSph) assuming all three potentials.
a All metal-rich globular clusters
b Metal-rich globular clusters excluding those associated with the Sagittarius dSph

1 and 2. In particular, old, metal-rich globular clus- cluster system with many accreted globular clusters, age
ters tend to have more negative specific orbital energies and metallicity should be correlated with each other but
and are therefore more tightly bound to the Milky Way inversely correlated with specific orbital energy. Our ob-
than young, metal-poor globular clusters. The metal- servation that old, metal-rich globular clusters are more
poor subsample prefers the full age–metallicity–specific tightly bound to the Milky Way than young, metal-poor
orbital energy relation over any of its nested submod- globular clusters confirms the importance of accretion
els. These results are independent of assumed poten- for the formation of the Milky Way’s globular cluster
tial, although the relationship is much stronger in the system. Our inference also supports previous propos-
McMillan17 potential from McMillan (2017) than either als that a large fraction of the Milky Way’s metal-poor
the default or scaled MWPotential2014 potentials from globular clusters were accreted (e.g., Forbes & Bridges
Bovy (2015). We find no clear age–metallicity–specific 2010; Forbes 2020; Kruijssen et al. 2019a,b, 2020; Mas-
orbital energy relation in the metal-rich subsample of 10 sari et al. 2019; Trujillo-Gomez et al. 2021).
globular clusters with [M/H] > −0.8 that are not asso- All but two of the metal rich clusters in our sample
ciated with the Sagittarius dSph galaxy. We also find have tightly bound orbits and old ages, broadly consis-
in this metal-rich subsample that an age–metallicity– tent with an in situ origin for metal-rich globular clus-
specific orbital energy relation is not clearly preferred ters. These two clusters—Terzan 7 and Pal 12—are both
over lower-dimensional models relating age, metallicity, young and on energetic orbits. They have been securely
or specific orbital energy. associated with the Sagittarius dSph galaxy though, and
As we argued in Section 1, accreted metal-rich globu- therefore have been accreted (e.g., Law & Majewski
lar clusters can only be formed in massive dSph galaxies 2010; Sohn et al. 2018). Because we do not observe a
that experience stronger dynamical friction and form strong age–metallicity–specific orbital energy for metal-
their first stars more quickly than lower-mass dSph rich globular clusters unassociated with the Sagittarius
galaxies. These two properties imply that in a globular dSph, our inferences do not support an accretion origin
 The Age–Metallicity–Specific Orbital Energy Relation for MW Globular Clusters 7

Table 4. Median Model Comparison F Statistics

Sample H0 : βτ = 0 H0 : βM = 0 H0 : βτ = βM = 0

F p F p F p
Default MWPotential2014
Complete 2.63 1.1 × 10−1 3.27 7.2 × 10−2 2.38 1.0 × 10−1
Metal-poor 5.71 2.1 × 10−2 4.66 3.7 × 10−2 3.58 3.7 × 10−2
Metal-richa 0.11 7.5 × 10−1 0.85 3.8 × 10−1 0.43 6.6 × 10−1
Metal-richb 3.19 1.1 × 10−1 0.85 3.9 × 10−1 1.62 2.6 × 10−1
Scaled MWPotential2014
Complete 2.17 1.5 × 10−1 4.40 4.1 × 10−2 2.65 8.0 × 10−2
Metal-poor 6.44 1.5 × 10−2 5.12 2.9 × 10−2 4.00 2.6 × 10−2
Metal-richa 0.00 9.8 × 10−1 1.27 2.9 × 10−1 0.68 5.3 × 10−1
Metal-richb 4.78 6.5 × 10−1 0.93 3.7 × 10−1 2.40 1.6 × 10−1
McMillan17
Complete 32.05 5.9 × 10−7 20.13 3.8 × 10−5 20.81 2.0 × 10−7
Metal-poor 20.89 4.2 × 10−5 16.38 2.2 × 10−4 12.90 4.3 × 10−5
Metal-richa 11.00 9.0 × 10−3 1.20 3.0 × 10−1 5.52 3.0 × 10−2
Metal-richb 0.55 4.8 × 10−1 0.09 7.7 × 10−1 0.49 6.3 × 10−1
Note—F statistics and associated p values comparing the full three parameter age–
metallicity–specific orbital energy relation derived for each potential to all possible lower-
dimensional submodels. The large F statistics and small p values indicate that the full
age–metallicity–specific orbital energy relation is superior to all lower-dimensional nested
models for the metal-poor subsample.
a All metal-rich globular clusters
b Metal-rich globular clusters excluding those associated with the Sagittarius dSph

for most metal-rich globular clusters. As there are only or theoretical model grid. This consistency supports the
10 metal-rich clusters unassociated with the Sagittarius robustness of our findings.
dSph in our sample, we are unable to choose between As a consequence of the simple physical explanation
gas-rich mergers or formation in gas-rich disks for the we outline to explain the age–metallicity–specific orbital
origin of the Milky Way’s metal-rich globular clusters. energy relation we observe in the Milky Way’s globular
In our analysis, we used relative ages from Marı́n- cluster system, we predict that the globular cluster sys-
Franch et al. (2009) calculated assuming the Carretta tems of other L∗ galaxies should evince similar relations.
& Gratton (1997) metallicity scale and the theoretical That is, we predict that other L∗ galaxies in isolation or
model grid from Dotter et al. (2007). To ensure that in low-density group environments like the Local Group
our results are insensitive to the assumed metallicity should have globular cluster systems in which old, metal-
scale or theoretical model grid underlying the Marı́n- rich globular clusters have more negative specific orbital
Franch et al. (2009) age inferences, we performed the energies and are therefore more tightly bound to their
same analysis using all of the relative ages from Marı́n- host halo than young, metal-poor globular clusters. We
Franch et al. (2009) calculated assuming both the Zinn argue that this is a natural outcome of galaxy formation
& West (1984) and Carretta & Gratton (1997) metal- in a ΛCDM universe.
licity scales and all of the Bertelli et al. (1994), Girardi While we infer the importance of accretion for the ori-
et al. (2000), Pietrinferni et al. (2004), and Dotter et al. gin of the Milky Way’s globular cluster system from an
(2007) theoretical model grids. In all cases, our analysis analysis of its ensemble properties, our approach can-
produces the same trends with comparable statistical not associate individual globular clusters with individ-
significance. This should not be surprising, as Marı́n- ual accretion events. Detailed simulations of Milky Way
Franch et al. (2009) found that their relative ages infer- analogs have shown how a globular cluster systems prop-
ences were insensitive to their choice of metallicity scale erties can be used to reveal that halo’s accretion history
8 Woody & Schlaufman

Metal Rich, Sgr dSph

Metal Rich, not Sgr dSph
Metal Poor

Specific Orbita
-4.0

-6.0

l Energy 10
-8.0

-10.0

h km 4
s2
1.2 -0.50

2
i
No 1.0
rm -1.00
aliz
ed
Ag h0.8 -1.50 /H]
e A ge
[M
12
.8
i 0.6 -2.00
Gy
r
Specific Orbita

-8.0

-10.0
l Energy 10

-12.0

-14.0
h

-16.0
km 4
s2

1.2 -0.50
2
i

No 1.0
rm -1.00
aliz
ed
Ag h0.8 -1.50 /H]
e A ge
[M
12
.8
i 0.6 -2.00
Gy
r
Specific Or

-5.0
bital Energ

-10.0

-15.0

-20.0
y 10
h

-25.0
4 km

1.2
s2

-0.50
2
i

No 1.0
rm -1.00
ali
ze
dA 0.8 H]
ge h -1.50 /
A
[M
12 ge
i 0.6
.8
Gy -2.00
r

Figure 1. Age–metallicity–specific orbital energy distribution of the Milky Way globular clusters listed in Tables 1 and 2.
We plot metal-poor clusters with [M/H] ≤ −0.8 as blue circles and metal-rich clusters with [M/H] > −0.8 as green triangles.
Downward-pointing triangles are metal-rich clusters associated with the Sagittarius dSph galaxy, while upward-pointing triangles
are metal-rich clusters unassociated with the Sagittarius dSph galaxy. Each plot is looking down the plane defined by the
median best-fit age–metallicity–specific orbital energy relation for the metal-poor subsample as reported in Table 3. Top:
specific orbital energies calculated assuming the default MWPotential2014 potential from Bovy (2015). Middle: specific orbital
energies calculated assuming the scaled MWPotential2014 potential. Bottom: specific orbital energies calculated assuming the
McMillan17 potential from McMillan (2017). It is visually clear in the McMillan17 potential that the metal-poor subsample lies
in a plane in age–metallicity–specific orbital energy space as indicated by the highly significant coefficients βτ and βM listed in
Table 3 and the small F test p values listed in Table 4.

(e.g., Kruijssen et al. 2019a,b, 2020; Pfeffer et al. 2020; associated with known accretion events. For example,
Trujillo-Gomez et al. 2021). When applied to the Milky Massari et al. (2019) employed this technique in combi-
Way, that formalism can quantify the properties of mas- nation with in-group age–metallicity relations to assign
sive satellite galaxies that have since merged with the candidate accreted globular clusters to specific progeni-
Milky Way (e.g., Myeong et al. 2018; Hughes et al. 2019; tors. However, our analysis shows that it is beneficial to
Kruijssen et al. 2019b, 2020; Pfeffer et al. 2020; Trujillo- use similar techniques on the entire cluster population,
Gomez et al. 2021). Moreover, efforts have been made thereby capturing the general imprint that accretion has
to match the kinematics of clusters to stellar streams
 The Age–Metallicity–Specific Orbital Energy Relation for MW Globular Clusters 9

1.2

1.1
i
12.8 Gyr

1.0
Age
h

0.9
Normalized Age

0.8
Metallicity Uncertainty
0.7
Metal Rich, Sgr dSph
0.6 Metal Rich, not Sgr dSph
Metal Poor
0.5
-2.0 -1.8 -1.6 -1.4 -1.2 -1.0 -0.8 -0.6 -0.4
[M/H]

-4.0
i
km2
s2
Specific Orbital Energy 104

-6.0
h

-8.0

-10.0
Metallicity Uncertainty

-12.0 Metal Rich, Sgr dSph

Metal Rich, not Sgr dSph
Metal Poor
-14.0
-2.0 -1.8 -1.6 -1.4 -1.2 -1.0 -0.8 -0.6 -0.4
[M/H]

-4.0
i
km2
s2
Specific Orbital Energy 104

-6.0
h

-8.0

-10.0
Median Normalized Age Uncertainty

-12.0 Metal Rich, Sgr dSph

Metal Rich, not Sgr dSph
Metal Poor
-14.0
0.6 0.7 0.8 0.9 1.0 1.1
h i
Age
Normalized Age 12.8 Gyr

Figure 2. Two-dimensional projections of the age–metallicity–specific orbital energy distribution of the 57 Milky Way globular
clusters listed in Tables 1 and 2 assuming the default MWPotential2014 potential (Bovy 2015). We plot metal-poor clusters
as blue circles and metal-rich clusters as green triangles. Downward-pointing triangles are metal-rich clusters associated with
the Sagittarius dSph galaxy, while upward-pointing triangles are metal-rich clusters unassociated with the Sagittarius dSph
galaxy. The sizes of the points corresponds to the relative values of the suppressed axis. In blue we indicate the two-dimensional
projections of the best-fit age–metallicity–specific orbital energy relation in the metal-poor subsample for each of the 100 Monte
Carlo trials. Top: metallicity–age relation. Middle: metallicity–specific orbital energy relation. Bottom: age–specific orbital
energy relation.
10 Woody & Schlaufman

1.2

1.1
i
12.8 Gyr

1.0
Age
h

0.9
Normalized Age

0.8
Metallicity Uncertainty
0.7
Metal Rich, Sgr dSph
0.6 Metal Rich, not Sgr dSph
Metal Poor
0.5
-2.0 -1.8 -1.6 -1.4 -1.2 -1.0 -0.8 -0.6 -0.4
[M/H]

-8.0
i
km2
s2

-10.0
Specific Orbital Energy 104
h

-12.0

-14.0

-16.0
Metallicity Uncertainty

-18.0
Metal Rich, Sgr dSph
-20.0 Metal Rich, not Sgr dSph
Metal Poor
-22.0
-2.0 -1.8 -1.6 -1.4 -1.2 -1.0 -0.8 -0.6 -0.4
[M/H]

-6.0
i

-8.0
km2
s2
Specific Orbital Energy 104

-10.0
h

-12.0

-14.0

-16.0
Median Normalized Age Uncertainty
-18.0
Metal Rich, Sgr dSph
-20.0 Metal Rich, not Sgr dSph
Metal Poor
-22.0
0.6 0.7 0.8 0.9 1.0 1.1
h i
Age
Normalized Age 12.8 Gyr

Figure 3. Two-dimensional projections of the age–metallicity–specific orbital energy distribution of the 57 Milky Way globular
clusters listed in Tables 1 and 2 assuming the scaled MWPotential2014 potential. We plot metal-poor clusters as blue circles
and metal-rich clusters as green triangles. Downward-pointing triangles are metal-rich clusters associated with the Sagittarius
dSph galaxy, while upward-pointing triangles are metal-rich clusters unassociated with the Sagittarius dSph galaxy. The sizes
of the points corresponds to the relative values of the suppressed axis. In blue we indicate the two-dimensional projections of
the best-fit age–metallicity–specific orbital energy relation in the metal-poor subsample for each of the 100 Monte Carlo trials.
Top: metallicity–age relation. Middle: metallicity–specific orbital energy relation. Bottom: age–specific orbital energy relation.
 The Age–Metallicity–Specific Orbital Energy Relation for MW Globular Clusters 11

1.2

1.1
i
12.8 Gyr

1.0
Age
h

0.9
Normalized Age

0.8
Metallicity Uncertainty
0.7
Metal Rich, Sgr dSph
0.6 Metal Rich, not Sgr dSph
Metal Poor
0.5
-2.0 -1.8 -1.6 -1.4 -1.2 -1.0 -0.8 -0.6 -0.4
[M/H]

-7.5
i

-10.0
km2
s2
Specific Orbital Energy 104

-12.5
h

-15.0

-17.5

Metallicity Uncertainty
-20.0

Metal Rich, Sgr dSph

-22.5
Metal Rich, not Sgr dSph
-25.0 Metal Poor

-2.0 -1.8 -1.6 -1.4 -1.2 -1.0 -0.8 -0.6 -0.4

[M/H]

-7.5
i

-10.0
km2
s2
Specific Orbital Energy 104

-12.5
h

-15.0

-17.5

-20.0 Median Normalized Age Uncertainty

Metal Rich, Sgr dSph

-22.5
Metal Rich, not Sgr dSph
-25.0 Metal Poor

0.6 0.7 0.8 0.9 1.0 1.1

h i
Age
Normalized Age 12.8 Gyr

Figure 4. Two-dimensional projections of the age–metallicity–specific orbital energy distribution of the 57 Milky Way globular
clusters listed in Tables 1 and 2 assuming McMillan17 potential (McMillan 2017). We plot metal-poor clusters as blue circles
and metal-rich clusters as green triangles. Downward-pointing triangles are metal-rich clusters associated with the Sagittarius
dSph galaxy, while upward-pointing triangles are metal-rich clusters unassociated with the Sagittarius dSph galaxy. The sizes
of the points corresponds to the relative values of the suppressed axis. In blue we indicate the two-dimensional projections of
the best-fit age–metallicity–specific orbital energy relation in the metal-poor subsample for each of the 100 Monte Carlo trials.
Top: metallicity–age relation. Middle: metallicity–specific orbital energy relation. Bottom: age–specific orbital energy relation.
 12 Woody & Schlaufman

Table 5. Median Model Comparison BIC Statistics they are in the McMillan17 potential, though this dif-
ference is diminished in the scaled MWPotential2014.
Sample Full Model βτ = 0 βM = 0 βτ = βM = 0
This results in a relatively larger spread in orbital en-
Default MWPotential2014 ergy, which might explain the weaker age–metallicity–
Complete 241 240 240 238 specific orbital energy relation. Nevertheless, the age–
Metal-poor 188 190 189 188
metallicity–specific orbital energy relation has the same
Metal-richa 55 53 54 52
form in all three potentials and we are confident its ex-
Metal-richb 43 44 41 42
istence is robust to the choice of potential.
Scaled MWPotential2014
Complete 247 245 248 244
Metal-poor 192 194 193 192
Metal-richa 57 54 56 53
Metal-richb 42 45 41 43
McMillan17
Complete 302 324 316 326
Metal-poor 236 251 247 250
Metal-richa 69 76 68 74
Metal-richb 55 54 53 52

Note—BIC values comparing the full three parameter age–metallicity–specific

orbital energy relation derived for each potential to all possible lower-
dimensional submodels. The small BIC values indicate that the full age–
metallicity–specific orbital energy relation is superior to all lower-dimensional
nested models for both the complete sample and its metal-poor subsample in
the McMillan17 potential.
a All metal-rich globular clusters 5. CONCLUSION
b Metal-rich globular clusters excluding those associated with the Sagittarius We find that in the metal-poor subsample of 45 Milky
dSph Way globular clusters with [M/H] ≤ −0.8 in Tables 1
and 2, relatively young or metal-poor globular clusters
are weakly bound to the Milky Way, while relatively
old or metal-rich globular clusters are tightly bound to
had on the Milky Way’s globular cluster system’s age– the Galaxy. We argue that this relationship is naturally
metallicity–specific orbital energy relation. explained by the accretion of globular clusters from now-
Because proper motions are more difficult to infer for disrupted dwarf galaxies. We propose that this obser-
distant clusters, our input sample lacks outer halo clus- vation is a consequence of the combined effects of the
ters. This bias is unlikely to affect our conclusion, as stellar mass–metallicity relationship for dwarf galaxies,
the only plausible explanation for weakly bound outer the dependence of dynamical friction on mass, and the
halo clusters is accretion from dwarf galaxies. Our in- relationship between a galaxy’s dynamical time and its
put sample also lacks highly extincted disk and bulge central density. Accreted metal-rich globular clusters
clusters. This bias is unlikely to affect our conclusion can only form in massive dSph galaxies that can at-
either, as disk/bulge clusters are metal-rich and would tain high metallicities and will be strongly affected by
not appear in our metal-poor subsample even if they dynamical friction. These more massive dSph galaxies
were present in our input sample. capable of forming metal-rich clusters will also produce
The age–metallicity–specific orbital energy relation we their oldest stellar populations before lower-mass dSph
observe in our metal-poor subsample of globular clusters galaxies that form in lower σ peaks in the universe’s pri-
is much stronger in the McMillan17 potential than in mordial matter density distribution. These lower-mass
the MWPotential2014 potential. As discussed in Section dSph galaxies will be much less affected by dynamical
3, the default MWPotential2014 potential represents a friction and can only produce low-metallicity globular
comparatively smaller, lower-mass Milky Way than the clusters. Because these properties are independent of
McMillan17 potential. In addition, the more massive the specific accretion history of the Milky Way, we as-
McMillan17 potential appears to be a better match to sert that they are a natural outcome of galaxy formation
Gaia DR2-informed Milky Way virial mass inferences. in a ΛCDM universe. We predict that the globular clus-
Under the default MWPotential2014 potential, most of ter systems of other L∗ galaxies in isolation or in low-
the weakly bound clusters are relatively less bound than density groups like the Local Group will lie in a plane
in age–metallicity–specific orbital energy space.
 The Age–Metallicity–Specific Orbital Energy Relation for MW Globular Clusters 13

ACKNOWLEDGMENTS
We thank Brendan Griffen and David Nataf for in-
sightful suggestions that improved our analyses. 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 Mul-
tilateral Agreement. This research has made use of
NASA’s Astrophysics Data System Bibliographic Ser-
vices. This research has made use of the SIMBAD
database, operated at CDS, Strasbourg, France (Wenger
et al. 2000). This research has made use of the VizieR
catalogue access tool, CDS, Strasbourg, France. The
original description of the VizieR service was published
in A&AS 143, 23 (Ochsenbein et al. 2000).

Software: galpy (Bovy 2015), numpy (van der Walt

et al. 2011), pandas (McKinney 2010), statsmodels
(Seabold & Perktold 2010)
14 Woody & Schlaufman

REFERENCES
Arenou, F., Luri, X., Babusiaux, C., et al. 2018, A&A, 616, Girardi, L., Bressan, A., Bertelli, G., & Chiosi, C. 2000,
A17 A&AS, 141, 371
Ashman, K. M., & Zepf, S. E. 1992, ApJ, 384, 50 Gravity Collaboration, Abuter, R., Amorim, A., et al. 2018,
Baumgardt, H., Hilker, M., Sollima, A., & Bellini, A. 2019, A&A, 615, L15
MNRAS, 482, 5138 Griffen, B. F., Drinkwater, M. J., Thomas, P. A., Helly,
Bertelli, G., Bressan, A., Chiosi, C., Fagotto, F., & Nasi, E. J. C., & Pimbblet, K. A. 2010, MNRAS, 405, 375
1994, A&AS, 106, 275 Hambly, N. C., Cropper, M., Boudreault, S., et al. 2018,
Binney, J., & Tremaine, S. 2008, Galactic Dynamics: A&A, 616, A15
Second Edition Harris, W. E. 1996, AJ, 112, 1487
Bjork, S. R., & Chaboyer, B. 2006, ApJ, 641, 1102 Hasselquist, S., Shetrone, M., Smith, V., et al. 2017, ApJ,
Bland-Hawthorn, J., & Gerhard, O. 2016, ARA&A, 54, 529 845, 162
Bovy, J. 2015, ApJS, 216, 29 Hughes, M. E., Pfeffer, J., Martig, M., et al. 2019, MNRAS,
Buonanno, R., Corsi, C. E., Castellani, M., et al. 1999, AJ, 482, 2795
118, 1671 Jurić, M., Ivezić, Ž., Brooks, A., et al. 2008, ApJ, 673, 864
Buonanno, R., Corsi, C. E., Zinn, R., et al. 1998, ApJL, Keller, B. W., Kruijssen, J. M. D., Pfeffer, J., et al. 2020,
501, L33 MNRAS, 495, 4248
Carretta, E., Bragaglia, A., Gratton, R. G., et al. 2014, Kim, J.-h., Ma, X., Grudić, M. Y., et al. 2018, MNRAS,
A&A, 561, A87 474, 4232
Carretta, E., & Gratton, R. G. 1997, A&AS, 121, 95 Kirby, E. N., Cohen, J. G., Guhathakurta, P., et al. 2013,
Carretta, E., Bragaglia, A., Gratton, R. G., et al. 2010, ApJ, 779, 102
A&A, 520, A95 Kravtsov, A. V., & Gnedin, O. Y. 2005, ApJ, 623, 650
Chaboyer, B., Fenton, W. H., Nelan, J. E., Patnaude, D. J., Kruijssen, J. M. D. 2015, MNRAS, 454, 1658
& Simon, F. E. 2001, ApJ, 562, 521 —. 2019, MNRAS, 486, L20
Choksi, N., Gnedin, O. Y., & Li, H. 2018, MNRAS, 480, Kruijssen, J. M. D., Pfeffer, J. L., Crain, R. A., & Bastian,
2343 N. 2019a, MNRAS, 486, 3134
Chou, M.-Y., Majewski, S. R., Cunha, K., et al. 2007, ApJ, Kruijssen, J. M. D., Pfeffer, J. L., Reina-Campos, M.,
670, 346 Crain, R. A., & Bastian, N. 2019b, MNRAS, 486, 3180
Cohen, J. G. 2004, AJ, 127, 1545 Kruijssen, J. M. D., Pfeffer, J. L., Chevance, M., et al.
Crowley, C., Kohley, R., Hambly, N. C., et al. 2016, A&A, 2020, MNRAS, 498, 2472
595, A6 Larsen, S. S., Brodie, J. P., & Strader, J. 2012, A&A, 546,
Dotter, A., Chaboyer, B., Jevremović, D., et al. 2007, AJ, A53
134, 376 Law, D. R., & Majewski, S. R. 2010, ApJ, 714, 229
Eadie, G., & Jurić, M. 2019, The Astrophysical Journal, Leaman, R., VandenBerg, D. A., & Mendel, J. T. 2013,
875, 159 MNRAS, 436, 122
Erkal, D., Belokurov, V., Laporte, C. F. P., et al. 2019, Letarte, B., Hill, V., Jablonka, P., et al. 2006, A&A, 453,
MNRAS, 487, 2685 547
Fabricius, C., Bastian, U., Portell, J., et al. 2016, A&A, Li, H., & Gnedin, O. Y. 2014, ApJ, 796, 10
595, A3 —. 2019, MNRAS, 486, 4030
Forbes, D. A. 2020, MNRAS, 493, 847 Lindegren, L., Hernández, J., Bombrun, A., et al. 2018,
Forbes, D. A., & Bridges, T. 2010, MNRAS, 404, 1203 A&A, 616, A2
Forbes, D. A., Bastian, N., Gieles, M., et al. 2018, Luri, X., Brown, A. G. A., Sarro, L. M., et al. 2018, A&A,
Proceedings of the Royal Society of London Series A, 616, A9
474, 20170616 Mackey, A. D., & Gilmore, G. F. 2004, MNRAS, 355, 504
Gaia Collaboration, Prusti, T., de Bruijne, J. H. J., et al. Marı́n-Franch, A., Aparicio, A., Piotto, G., et al. 2009,
2016, A&A, 595, A1 ApJ, 694, 1498
Gaia Collaboration, Helmi, A., van Leeuwen, F., et al. Massari, D., Koppelman, H. H., & Helmi, A. 2019, A&A,
2018a, A&A, 616, A12 630, L4
Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. McKinney, W. 2010, in Proceedings of the 9th Python in
2018b, A&A, 616, A1 Science Conference (SciPy)
 The Age–Metallicity–Specific Orbital Energy Relation for MW Globular Clusters 15

McMillan, P. J. 2017, MNRAS, 465, 76 Rutledge, G. A., Hesser, J. E., Stetson, P. B., et al. 1997b,
Miyamoto, M., & Nagai, R. 1975, PASJ, 27, 533 PASP, 109, 883
Mo, H., van den Bosch, F. C., & White, S. 2010, Galaxy Sbordone, L., Bonifacio, P., Buonanno, R., et al. 2007,
Formation and Evolution A&A, 465, 815
Mottini, M., Wallerstein, G., & McWilliam, A. 2008, AJ, Schönrich, R., Binney, J., & Dehnen, W. 2010, MNRAS,
136, 614 403, 1829
Seabold, S., & Perktold, J. 2010, in 9th Python in Science
Muratov, A. L., & Gnedin, O. Y. 2010, ApJ, 718, 1266
Conference
Myeong, G. C., Evans, N. W., Belokurov, V., Sanders,
Searle, L., & Zinn, R. 1978, ApJ, 225, 357
J. L., & Koposov, S. E. 2018, ApJL, 863, L28
Simon, J. D. 2019, ARA&A, 57, 375
Navarro, J. F., Frenk, C. S., & White, S. D. M. 1996, ApJ,
Sohn, S. T., Watkins, L. L., Fardal, M. A., et al. 2018, ApJ,
462, 563
862, 52
Ochsenbein, F., Bauer, P., & Marcout, J. 2000, A&AS, 143,
Trujillo-Gomez, S., Kruijssen, J. M. D., Reina-Campos, M.,
23 et al. 2021, MNRAS, 503, 31
Pascale, R., Posti, L., Nipoti, C., & Binney, J. 2018, van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011,
MNRAS, 480, 927 Computing in Science and Engineering, 13, 22
Patel, S. G., van Dokkum, P. G., Franx, M., et al. 2013a, VandenBerg, D. A., Brogaard, K., Leaman, R., &
ApJ, 766, 15 Casagrand e, L. 2013, ApJ, 775, 134
Patel, S. G., Fumagalli, M., Franx, M., et al. 2013b, ApJ, Vasiliev, E., & Belokurov, V. 2020, MNRAS, 497, 4162
778, 115 Wagner-Kaiser, R., Mackey, D., Sarajedini, A., et al. 2017,
Pfeffer, J., Kruijssen, J. M. D., Crain, R. A., & Bastian, N. MNRAS, 471, 3347
2018, MNRAS, 475, 4309 Wang, W., Han, J., Cautun, M., Li, Z., & Ishigaki, M. N.
Pfeffer, J. L., Trujillo-Gomez, S., Kruijssen, J. M. D., et al. 2020, Science China Physics, Mechanics, and Astronomy,
2020, MNRAS, 499, 4863 63, 109801
Pietrinferni, A., Cassisi, S., Salaris, M., & Castelli, F. 2004, Wenger, M., Ochsenbein, F., Egret, D., et al. 2000, A&AS,
ApJ, 612, 168 143, 9
Renaud, F., Agertz, O., & Gieles, M. 2017, MNRAS, 465, West, M. J., Côté, P., Marzke, R. O., & Jordán, A. 2004,
3622 Nature, 427, 31
Wright, E. L. 2006, PASP, 118, 1711
Rutledge, G. A., Hesser, J. E., & Stetson, P. B. 1997a,
Zinn, R., & West, M. J. 1984, ApJS, 55, 45
PASP, 109, 907
 16

Table 1. Input Astrometric, Distance, and Radial Velocity Data

NGC Alternate α δ µα µδ $ σαδ σα$ σδ$ d vr

Name Name [deg] [deg] [mas yr−1 ] [mas yr−1 ] [mas] [kpc] [km s−1 ]

NGC 104 47 Tuc 6.0194 −72.0821 5.2477+0.0016

−0.0016 −2.5189+0.0015
−0.0015 0.1959+0.0002
−0.0002 −0.06 −0.01 −0.01 4.3+0.1
−0.1 −18.7+0.2
−0.2
NGC 288 ... 13.1879 −26.5858 4.2385+0.0035
−0.0035 −5.6470+0.0026
−0.0026 0.1401+0.0021
−0.0021 0.25 0.15 −0.13 8.1+0.2
−0.2 −46.6+0.4
−0.4
NGC 362 ... 15.8099 −70.8489 6.6954+0.0045
−0.0045 −2.5184+0.0034
−0.0034 0.0788+0.0012
−0.0012 −0.09 −0.04 −0.12 8.3+0.2
−0.2 223.5+0.5
−0.5
NGC 1261 ... 48.067543 −55.216224 1.6900+0.0400
−0.0400 −2.1100+0.0400
−0.0400 ... ... ... ... 16.0+0.4
−0.4 68.2+4.6
−4.6
NGC 1851 ... 78.528 −40.0456 2.1308+0.0037
−0.0037 −0.6220+0.0040
−0.0040 0.0298+0.0011
−0.0011 −0.09 0.06 −0.07 12.2+0.3
−0.3 320.9+1.0
−1.0
NGC 2298 ... 102.2464 −36.0046 3.2762+0.0060
−0.0060 −2.1913+0.0061
−0.0061 0.0791+0.0019
−0.0019 0.07 0.08 −0.07 10.6+0.2
−0.2 148.9+1.2
−1.2
NGC 2808 ... 138.0071 −64.8645 1.0032+0.0032
−0.0032 0.2785+0.0032
−0.0032 0.0560+0.0006
−0.0006 −0.08 0.05 −0.01 9.3+0.2
−0.2 93.6+2.4
−2.4
NGC 3201 ... 154.3987 −46.4125 8.3344+0.0021
−0.0021 −1.9895+0.0020
−0.0020 0.1724+0.0006
−0.0006 0.12 0.04 −0.02 5.1+0.1
−0.1 494.0+0.2
−0.2
NGC 4147 ... 182.52626 18.542638 −1.7800+0.0400
−0.0400 −2.1000+0.0400
−0.0400 ... ... ... ... 18.8+0.4
−0.4 183.2+0.7
−0.7
NGC 4590 M 68 189.8651 −26.7454 −2.7640+0.0050
−0.0050 1.7916+0.0039
−0.0039 0.0664+0.0025
−0.0025 −0.29 −0.0 0.13 10.1+0.2
−0.2 −95.2+0.4
−0.4
NGC 4833 ... 194.8978 −70.8718 −8.3147+0.0036
−0.0036 −0.9366+0.0029
−0.0029 0.1163+0.0010
−0.0010 0.06 0.05 0.11 5.9+0.1
−0.1 200.2+1.2
−1.2
NGC 5024 M 53 198.2262 18.1661 −0.1466+0.0045
−0.0045 −1.3514+0.0032
−0.0032 0.0143+0.0018
−0.0018 −0.28 −0.12 0.08 18.4+0.4
−0.4 −79.1+4.1
−4.1
NGC 5053 ... 199.1124 17.7008 −0.3591+0.0071
−0.0071 −1.2586+0.0048
−0.0048 0.0064+0.0040
−0.0040 −0.32 0.13 −0.17 16.2+0.4
−0.4 44.0+0.4
−0.4
NGC 5139 ω Cen 201.7876 −47.4515 −3.1925+0.0022
−0.0022 −6.7445+0.0019
−0.0019 0.1237+0.0011
−0.0011 −0.03 −0.04 0.17 5.1+0.1
−0.1 232.3+0.5
−0.5
NGC 5272 M3 205.5486 28.376 −0.1127+0.0029
−0.0029 −2.6274+0.0022
−0.0022 0.0265+0.0010
−0.0010 −0.03 −0.01 −0.05 10.0+0.2
−0.2 −148.6+0.4
−0.4
NGC 5286 ... 206.6136 −51.3723 0.1836+0.0076
−0.0076 −0.1477+0.0068
−0.0068 0.0168+0.0025
−0.0025 −0.01 −0.02 0.08 10.7+0.2
−0.2 58.3+3.2
−3.2
NGC 5466 ... 211.3614 28.5331 −5.4044+0.0042
−0.0042 −0.7907+0.0041
−0.0041 0.0210+0.0021
−0.0021 0.07 0.04 0.15 16.6+0.4
−0.4 107.7+0.3
−0.3
NGC 5904 M5 229.6394 2.0766 4.0613+0.0032
−0.0032 −9.8610+0.0029
−0.0029 0.1135+0.0010
−0.0010 0.03 −0.07 0.09 7.3+0.2
−0.2 52.1+0.5
−0.5
NGC 5927 ... 232.0065 −50.6694 −5.0470+0.0060
−0.0060 −3.2325+0.0055
−0.0055 0.0996+0.0021
−0.0021 −0.08 −0.0 −0.01 7.4+0.2
−0.2 −115.7+3.1
−3.1
NGC 5986 ... 236.5211 −37.7826 −4.2217+0.0084
−0.0084 −4.5515+0.0065
−0.0065 0.0718+0.0031
−0.0031 −0.18 −0.1 −0.01 10.3+0.2
−0.2 88.9+3.7
−3.7
Woody & Schlaufman

NGC 6093 M 80 244.2564 −22.9723 −2.9469+0.0090

−0.0090 −5.5613+0.0073
−0.0073 0.0558+0.0030
−0.0030 0.01 −0.1 0.06 8.7+0.2
−0.2 9.3+3.1
−3.1
NGC 6121 M4 245.8976 −26.5279 −12.4956+0.0033
−0.0033 −18.9789+0.0030
−0.0030 0.5001+0.0007
−0.0007 0.08 −0.08 0.03 2.2+0.1
−0.1 70.2+0.3
−0.3
NGC 6101 ... 246.450485 −72.202194 1.7300+0.0500
−0.0500 −0.4700+0.0500
−0.0500 ... ... ... ... 15.1+0.3
−0.4 361.4+1.7
−1.7
NGC 6144 ... 246.8061 −26.0301 −1.7646+0.0085
−0.0085 −2.6371+0.0063
−0.0063 0.0668+0.0040
−0.0040 0.08 −0.15 0.15 10.1+0.2
−0.2 188.9+1.1
−1.1
NGC 6171 M 107 248.135 −13.057 −1.9359+0.0064
−0.0064 −5.9487+0.0048
−0.0048 0.1480+0.0026
−0.0026 0.03 −0.13 0.17 6.3+0.1
−0.1 −33.8+0.3
−0.3
NGC 6205 M 13 250.4217 36.4596 −3.1762+0.0027
−0.0027 −2.5876+0.0030
−0.0030 0.0801+0.0007
−0.0007 0.19 −0.04 0.04 7.0+0.2
−0.2 −246.6+0.9
−0.9
NGC 6218 M 12 251.8101 −1.951 −0.1577+0.0040
−0.0040 −6.7683+0.0027
−0.0027 0.1563+0.0013
−0.0013 0.29 −0.06 0.11 4.7+0.1
−0.1 −42.1+0.6
−0.6
NGC 6254 M 10 254.2861 −4.0981 −4.7031+0.0039
−0.0039 −6.5285+0.0027
−0.0027 0.1511+0.0014
−0.0014 0.21 −0.07 0.1 4.3+0.1
−0.1 75.8+1.0
−1.0
NGC 6304 ... 258.637 −29.4816 −3.9478+0.0095
−0.0095 −1.1248+0.0069
−0.0069 0.1077+0.0034
−0.0034 0.16 −0.04 0.13 6.0+0.1
−0.1 −107.3+3.6
−3.6
NGC 6341 M 92 259.2821 43.1352 −4.9367+0.0040
−0.0040 −0.5559+0.0040
−0.0040 0.0564+0.0008
−0.0008 0.11 −0.04 −0.0 8.1+0.2
−0.2 −121.6+1.5
−1.5
NGC 6352 ... 261.3739 −48.427 −2.1889+0.0046
−0.0046 −4.4209+0.0036
−0.0036 0.1543+0.0018
−0.0018 0.22 −0.12 0.13 5.6+0.1
−0.1 −120.9+3.0
−3.0
NGC 6366 ... 261.9393 −5.0752 −0.3835+0.0054
−0.0054 −5.1309+0.0044
−0.0044 0.2292+0.0022
−0.0022 0.27 −0.08 0.14 3.6+0.1
−0.1 −122.3+0.5
−0.5
NGC 6362 ... 262.9772 −67.0492 −5.5014+0.0028
−0.0028 −4.7417+0.0032
−0.0032 0.0974+0.0011
−0.0011 0.06 −0.06 0.06 7.5+0.2
−0.2 −13.1+0.6
−0.6
NGC 6388 ... 264.0654 −44.7423 −1.3548+0.0072
−0.0072 −2.7144+0.0061
−0.0061 0.0482+0.0034
−0.0034 0.11 −0.12 0.18 11.5+0.3
−0.3 81.2+1.2
−1.2
NGC 6397 ... 265.1697 −53.6773 3.2908+0.0026
−0.0026 −17.5908+0.0025
−0.0025 0.3781+0.0007
−0.0007 0.1 −0.05 0.07 2.2+0.1
−0.1 18.9+0.1
−0.1

Table 1 continued
 Table 1 (continued)

NGC Alternate α δ µα µδ $ σαδ σα$ σδ$ d vr

Name Name [deg] [deg] [mas yr−1 ] [mas yr−1 ] [mas] [kpc] [km s−1 ]

NGC 6441 ... 267.5540 −37.066 −2.5394+0.0070

−0.0070 −5.3010+0.0057
−0.0057 0.0403+0.0037
−0.0037 0.17 −0.18 0.17 9.7+0.2
−0.2 18.3+1.7
−1.7
NGC 6496 ... 269.7677 −44.2660 −3.0290+0.0057
−0.0057 −9.1971+0.0050
−0.0050 0.0803+0.0031
−0.0031 0.07 −0.15 0.17 11.6+0.3
−0.3 −112.7+5.7
−5.7
NGC 6535 ... 270.959 −0.2953 −4.2101+0.0115
−0.0115 −2.9461+0.0108
−0.0108 0.1294+0.0047
−0.0047 0.1 0.1 −0.02 6.8+0.2
−0.2 −215.1+0.5
−0.5
NGC 6541 ... 271.9827 −43.7144 0.2762+0.0054
−0.0054 −8.7659+0.0048
−0.0048 0.1139+0.0025
−0.0025 −0.03 −0.04 0.08 7.4+0.2
−0.2 −156.2+2.7
−2.7
NGC 6584 ... 274.656646 −52.215778 −0.0500+0.0100
−0.0100 −7.2200+0.0100
−0.0100 ... −0.15 ... ... 13.0+0.3
−0.3 222.9+15.0
−15.0
NGC 6637 M 69 277.8342 −32.3565 −5.0669+0.0104
−0.0104 −5.8017+0.0094
−0.0094 0.0746+0.0032
−0.0032 0.23 −0.01 −0.05 8.2+0.2
−0.2 39.9+2.8
−2.8
NGC 6656 M 22 279.1048 −23.9102 9.8019+0.0036
−0.0036 −5.5643+0.0034
−0.0034 0.2602+0.0009
−0.0009 0.15 −0.05 0.02 3.2+0.1
−0.1 −148.9+0.4
−0.4
NGC 6681 M 70 280.802 −32.2892 1.3853+0.0076
−0.0076 −4.7174+0.0065
−0.0065 0.1096+0.0038
−0.0038 0.26 −0.24 0.1 8.7+0.2
−0.2 218.6+1.2
−1.2
NGC 6723 ... 284.888123 −36.632248 1.0000+0.0100
−0.0100 −2.4200+0.0100
−0.0100 ... 0.11 ... ... 8.6+0.2
−0.2 −94.5+3.6
−3.6
NGC 6752 ... 287.7175 −59.9833 −3.1908+0.0018
−0.0018 −4.0347+0.0020
−0.0020 0.2310+0.0011
−0.0011 0.19 −0.29 0.03 3.9+0.1
−0.1 −24.5+1.9
−1.9
NGC 6779 M 56 289.148 30.184 −2.0092+0.0051
−0.0051 1.6553+0.0056
−0.0056 0.0702+0.0015
−0.0015 0.03 −0.05 −0.03 9.9+0.2
−0.2 −135.7+0.8
−0.8
Terzan 7 ... 289.432983 −34.657722 −3.0400+0.0600
−0.0600 −1.7100+0.0500
−0.0500 ... ... ... ... 23.0+0.5
−0.5 166.0+4.0
−4.0
Arp 2 ... 292.183807 −30.355638 −2.4000+0.0400
−0.0400 −1.5400+0.0300
−0.0300 ... ... ... ... 27.6+0.6
−0.6 115.0+10.0
−10.0
NGC 6809 M 55 295.0046 −30.9621 −3.4017+0.0031
−0.0031 −9.2642+0.0028
−0.0028 0.1707+0.0011
−0.0011 0.18 −0.03 0.0 5.3+0.1
−0.1 174.8+0.4
−0.4
Terzan 8 ... 295.435028 −33.999474 −2.9100+0.0800
−0.0800 −1.6300+0.0600
−0.0600 ... ... ... ... 25.4+0.6
−0.6 130.0+8.0
−8.0
NGC 6838 M 71 298.4427 18.779 −3.3842+0.0027
−0.0027 −2.6528+0.0028
−0.0028 0.2252+0.0010
−0.0010 0.11 −0.11 −0.01 3.8+0.1
−0.1 −22.9+0.2
−0.2
NGC 6934 ... 308.547393 7.404472 −2.6700+0.0400
−0.0400 −4.5200+0.0500
−0.0500 ... ... ... ... 15.2+0.3
−0.4 −411.4+1.6
−1.6
NGC 6981 M 72 313.3662 −12.5386 −1.2488+0.0089
−0.0089 −3.3117+0.0068
−0.0068 0.0225+0.0063
−0.0063 0.26 −0.38 −0.13 16.8+0.4
−0.4 −345.1+3.7
−3.7
NGC 7078 M 15 322.4949 12.1661 −0.6238+0.0041
−0.0041 −3.7960+0.0039
−0.0039 0.0568+0.0014
−0.0014 −0.04 −0.02 −0.15 10.2+0.2
−0.2 −107.5+0.3
−0.3
NGC 7089 M2 323.3497 −0.8177 3.4911+0.0077
−0.0077 −2.1501+0.0071
−0.0071 0.0591+0.0035
−0.0035 −0.04 −0.14 −0.14 11.4+0.3
−0.3 −5.3+2.0
−2.0
NGC 7099 M 30 325.0888 −23.1792 −0.7017+0.0063
−0.0063 −7.2218+0.0055
−0.0055 0.0746+0.0040
−0.0040 0.3 −0.29 −0.27 7.9+0.2
−0.2 −184.2+1.0
−1.0
Pal 12 ... 326.661804 −21.252611 −3.0600+0.0500
−0.0500 −3.3200+0.0500
−0.0500 ... ... ... ... 18.7+0.4
−0.4 27.8+1.5
−1.5

Note—Input globular cluster astrometry, distance, and radial velocity data. While most cluster astrometric data come from Gaia Collaboration
et al. (2018a), the data for NGC 6584 and NGC 6723 come from Baumgardt et al. (2019) while the data for NGC 1261, NGC 4147, NGC 6101,
Terzan 7, Arp 2, Terzan 8, NGC 6934, and Pal 12 come from Sohn et al. (2018). Distance and radial velocity data come from the December 2010
revision of the Harris (1996) compilation.
The Age–Metallicity–Specific Orbital Energy Relation for MW Globular Clusters
17
 18

Table 2. Age, Metallicity, and Specific Orbital Energy Inferences

NGC Alternate τ [M/H] SOEMW14 SOEScaledMW14 SOEMc17

2 2
Name Name [104 km s−2 ] [104 km s−2 ] [104 km2 s−2 ]
NGC 104 47 Tuc 1.02+0.07
−0.07 −0.64+0.05
−0.05 −7.64+0.06
−0.05 −11.82+0.05
−0.05 −16.62+0.05
−0.05
NGC 288 ... 0.83+0.03
−0.03 −0.92+0.05
−0.05 −9.04+0.08
−0.04 −13.06+0.07
−0.05 −16.00+0.11
−0.07
NGC 362 ... 0.81+0.04
−0.04 −0.87+0.05
−0.05 −8.66+0.16
−0.18 −12.84+0.17
−0.17 −16.20+0.17
−0.20
NGC 1261 ... 0.80+0.04
−0.04 −0.86+0.05
−0.05 −8.20+0.07
−0.06 −12.04+0.06
−0.08 −13.43+0.13
−0.09
NGC 1851 ... 0.78+0.04
−0.04 −0.81+0.05
−0.05 −7.57+0.15
−0.10 −11.34+0.14
−0.09 −13.47+0.14
−0.09
NGC 2298 ... 0.99+0.05
−0.05 −1.49+0.05
−0.05 −8.04+0.19
−0.10 −11.85+0.17
−0.13 −14.15+0.22
−0.15
NGC 2808 ... 0.85+0.03
−0.03 −0.89+0.05
−0.05 −8.12+0.08
−0.07 −12.24+0.13
−0.08 −15.33+0.11
−0.10
NGC 3201 ... 0.80+0.03
−0.03 −1.02+0.05
−0.05 −3.13+0.19
−0.25 −7.21+0.19
−0.25 −11.01+0.23
−0.26
NGC 4147 ... 0.89+0.04
−0.04 −1.28+0.05
−0.05 −7.34+0.04
−0.04 −11.08+0.06
−0.04 −12.33+0.11
−0.09
NGC 4590 M 68 0.90+0.04
−0.04 −1.78+0.05
−0.05 −3.94+0.56
−0.59 −8.17+0.71
−0.56 −11.41+0.51
−0.55
NGC 4833 ... 0.98+0.05
−0.05 −1.49+0.05
−0.05 −9.37+0.10
−0.09 −13.58+0.08
−0.08 −18.66+0.12
−0.11
NGC 5024 M 53 0.99+0.05
−0.05 −1.64+0.05
−0.05 −7.44+0.88
−0.53 −11.04+0.97
−0.70 −12.26+0.90
−0.53
NGC 5053 ... 0.96+0.04
−0.04 −1.76+0.05
−0.05 −7.75+0.73
−0.60 −11.68+0.92
−0.60 −12.97+0.75
−0.60
NGC 5139 ω Cen 0.90+0.05 −1.13+0.05 −9.29+0.11 −13.51+0.07 −18.60+0.12
Woody & Schlaufman

−0.05 −0.05 −0.06 −0.08 −0.07

NGC 5272 M3 0.89+0.04
−0.04 −1.12+0.05
−0.05 −7.55+0.13
−0.12 −11.61+0.13
−0.12 −14.06+0.13
−0.12
NGC 5286 ... 0.98+0.04
−0.04 −1.19+0.05
−0.05 −8.06+0.07
−0.08 −12.48+0.11
−0.07 −15.88+0.16
−0.13
NGC 5466 ... 1.06+0.05
−0.05 −1.98+0.05
−0.05 −3.35+0.65
−0.52 −7.27+0.44
−0.53 −9.00+0.73
−0.46
NGC 5904 M5 0.83+0.03
−0.03 −0.90+0.05
−0.05 −4.44+0.59
−0.39 −8.76+0.50
−0.41 −12.72+0.68
−0.43
NGC 5927 ... 0.99+0.07
−0.07 −0.50+0.05
−0.05 −7.18+0.09
−0.10 −11.62+0.10
−0.09 −18.37+0.10
−0.07
NGC 5986 ... 0.95+0.04
−0.04 −1.13+0.05
−0.05 −10.81+0.19
−0.16 −15.47+0.21
−0.12 −20.28+0.29
−0.22
NGC 6093 M 80 0.98+0.04
−0.04 −1.25+0.05
−0.05 −10.84+0.07
−0.05 −15.42+0.08
−0.07 −21.78+0.17
−0.12
NGC 6121 M4 0.98+0.05
−0.05 −0.83+0.05
−0.05 −9.72+0.02
−0.01 −13.93+0.02
−0.01 −19.69+0.05
−0.04
NGC 6101 ... 0.98+0.04
−0.04 −1.54+0.05
−0.05 −5.04+0.13
−0.11 −9.60+0.13
−0.11 −10.27+0.22
−0.26
NGC 6144 ... 1.08+0.05
−0.05 −1.52+0.05
−0.05 −9.58+0.18
−0.18 −14.31+0.15
−0.13 −19.82+0.18
−0.20

Table 2 continued
 Table 2 (continued)

NGC Alternate τ [M/H] SOEMW14 SOEScaledMW14 SOEMc17

Name Name [104 km2 s−2 ] [104 km2 s−2 ] [104 km2 s−2 ]
NGC 6171 M 107 1.09+0.06
−0.06 −0.81+0.05
−0.05 −10.00+0.09
−0.09 −14.44+0.08
−0.07 −21.27+0.11
−0.08
NGC 6205 M 13 0.91+0.04
−0.04 −1.11+0.05
−0.05 −9.32+0.10
−0.07 −13.51+0.10
−0.05 −17.34+0.11
−0.08
NGC 6218 M 12 0.99+0.03
−0.03 −0.92+0.05
−0.05 −9.01+0.15
−0.18 −13.32+0.18
−0.19 −19.47+0.17
−0.20
NGC 6254 M 10 0.89+0.04
−0.04 −1.03+0.05
−0.05 −8.82+0.10
−0.08 −13.14+0.11
−0.07 −19.37+0.12
−0.08
NGC 6304 ... 1.06+0.08
−0.08 −0.52+0.05
−0.05 −8.08+0.23
−0.21 −12.52+0.22
−0.22 −22.01+0.27
−0.34
NGC 6341 M 92 1.03+0.04
−0.04 −1.94+0.05
−0.05 −9.00+0.10
−0.11 −13.16+0.13
−0.10 −16.63+0.09
−0.11
NGC 6352 ... 0.99+0.07
−0.07 −0.56+0.05
−0.05 −7.69+0.17
−0.15 −12.06+0.13
−0.19 −20.12+0.18
−0.16
NGC 6366 ... 1.04+0.13
−0.13 −0.59+0.05
−0.05 −8.53+0.19
−0.10 −12.78+0.15
−0.11 −19.31+0.20
−0.10
NGC 6362 ... 1.06+0.05
−0.05 −0.85+0.05
−0.05 −8.97+0.08
−0.07 −13.36+0.07
−0.09 −19.02+0.09
−0.11
NGC 6388 ... 0.94+0.08
−0.08 −0.63+0.05
−0.05 −11.34+0.21
−0.16 −16.13+0.27
−0.18 −20.99+0.26
−0.29
NGC 6397 ... 0.99+0.04
−0.04 −1.54+0.05
−0.05 −8.24+0.04
−0.06 −12.44+0.05
−0.06 −18.30+0.04
−0.06
NGC 6441 ... 0.88+0.07
−0.07 −0.46+0.05
−0.05 −11.26+0.09
−0.08 −16.03+0.11
−0.08 −24.53+0.40
−0.33
NGC 6496 ... 0.97+0.07
−0.07 −0.56+0.05
−0.05 −5.86+0.60
−0.60 −10.94+0.69
−0.43 −15.40+0.77
−0.79
NGC 6535 ... 0.82+0.09
−0.09 −1.29+0.05
−0.05 −9.38+0.09
−0.09 −13.84+0.12
−0.11 −20.81+0.08
−0.08
NGC 6541 ... 1.01+0.04
−0.04 −1.31+0.05
−0.05 −7.52+0.11
−0.13 −12.05+0.15
−0.14 −20.74+0.12
−0.13
NGC 6584 ... 0.88+0.03
−0.03 −1.10+0.05
−0.05 −4.81+0.78
−0.75 −9.54+0.64
−0.73 −11.81+0.92
−0.91
NGC 6637 M 69 1.02+0.07
−0.07 −0.64+0.05
−0.05 −10.39+0.17
−0.17 −14.94+0.17
−0.16 −24.16+0.21
−0.17
NGC 6656 M 22 0.99+0.05
−0.05 −1.27+0.05
−0.05 −5.39+0.15
−0.16 −9.63+0.19
−0.19 −16.48+0.15
−0.17
NGC 6681 M 70 1.00+0.04
−0.04 −1.13+0.05
−0.05 −7.31+0.42
−0.35 −11.96+0.45
−0.34 −19.70+0.45
−0.34
NGC 6723 ... 1.02+0.05
−0.05 −0.82+0.05
−0.05 −9.40+0.43
−0.27 −14.01+0.42
−0.29 −21.21+0.49
−0.35
NGC 6752 ... 0.92+0.04
−0.04 −1.02+0.05
−0.05 −8.18+0.13
−0.09 −12.41+0.08
−0.10 −18.53+0.13
−0.10
NGC 6779 M 56 1.07+0.05
−0.05 −1.50+0.05
−0.05 −8.00+0.17
−0.14 −12.27+0.17
−0.16 −15.73+0.18
−0.17
Terzan 7 ... 0.57+0.04
−0.04 −0.42+0.05
−0.05 −9.39+0.21
−0.19 −14.83+0.22
−0.18 −9.46+0.34
−0.35
Arp 2 ... 0.85+0.06
−0.06 −1.23+0.05
−0.05 −9.59+0.43
−0.33 −14.87+0.35
−0.33 −9.10+0.44
−0.36
NGC 6809 M 55 0.96+0.04
−0.04 −1.32+0.05
−0.05 −8.19+0.05
−0.03 −12.56+0.05
−0.04 −18.88+0.05
−0.04
Terzan 8 ... 0.95+0.04
−0.04 −1.58+0.05
−0.05 −8.12+0.40
−0.37 −13.22+0.34
−0.39 −8.35+0.56
−0.35
The Age–Metallicity–Specific Orbital Energy Relation for MW Globular Clusters

NGC 6838 M 71 1.07+0.08

−0.08 −0.59+0.05
−0.05 −7.69+0.04
−0.05 −11.88+0.05
−0.05 −17.25+0.04
−0.04
19

Table 2 continued
 Table 2 (continued)
20

NGC Alternate τ [M/H] SOEMW14 SOEScaledMW14 SOEMc17

Name Name [104 km2 s−2 ] [104 km2 s−2 ] [104 km2 s−2 ]
NGC 6934 ... 0.87+0.04
−0.04 −1.10+0.05
−0.05 −4.97+0.20
−0.24 −9.37+0.22
−0.23 −10.43+0.33
−0.37
NGC 6981 M 72 0.85+0.02
−0.02 −0.99+0.05
−0.05 −7.58+0.67
−0.45 −12.04+0.71
−0.45 −12.66+0.79
−0.51
NGC 7078 M 15 1.01+0.04
−0.04 −1.80+0.05
−0.05 −9.24+0.09
−0.06 −13.43+0.09
−0.06 −16.15+0.13
−0.13
NGC 7089 M2 0.92+0.04
−0.04 −1.09+0.05
−0.05 −7.14+0.42
−0.39 −11.54+0.63
−0.42 −13.97+0.38
−0.39
NGC 7099 M 30 1.01+0.05
−0.05 −1.70+0.05
−0.05 −9.24+0.21
−0.12 −13.56+0.23
−0.13 −17.53+0.29
−0.17
Pal 12 ... 0.69+0.09
−0.09 −0.69+0.05
−0.05 −5.12+0.24
−0.26 −9.48+0.27
−0.24 −9.18+0.34
−0.35

Note—Table 2 is ordered by R.A. The first column is normalized age τ = age/12.8 Gyr from Marı́n-Franch
et al. (2009) calculated using the Dartmouth Stellar Evolution Program (Chaboyer et al. 2001; Bjork &
Chaboyer 2006; Dotter et al. 2007) assuming the metallicities presented in Rutledge et al. (1997a,b) on the
Carretta & Gratton (1997) scale. The second column metallicity [M/H] is from Rutledge et al. (1997a,b)
on the Carretta & Gratton (1997) scale. The third, fourth, and fifth columns are globular cluster specific
orbital energies calculated using the data in Table 1 assuming the default, MWPotential2014, the scaled
MWPotential2014 (Bovy 2015), and the McMillan17 (McMillan 2017) potentials respectively.
Woody & Schlaufman