THE MASS PROFILE OF THE GALAXY TO 80 kpc

Measuring the mass of the Galaxy is an astonishingly difficult task. Convenient tracers of disk rotation—stars and gas clouds—extend only to 20 kpc (e.g., Sofue et al. 2009). At larger radii, statistical analysis of radial velocities must be used. Traditional tracers at distances up to 100 kpc include globular clusters and dwarf satellite galaxies (e.g., Kochanek 1996; Wilkinson & Evans 1999; Sakamoto et al. 2003). Recently, Battaglia et al. (2005) derived the radial velocity dispersion profile to 120 kpc using a combined sample of globular clusters, satellite galaxies, and halo red giant stars. They found the dispersion falling from ∼120 km s⁻¹ to ∼50 km s⁻¹ between 20 and 120 kpc. In contrast, Xue et al. (2008) assembled a large sample of blue horizontal-branch (BHB) stars from the Sloan Digital Sky Survey Data Release 6 (SDSS DR6) and found a much flatter profile between 20 and 60 kpc. Here, we use a new spectroscopic survey (Brown et al. 2010) aimed at finding hypervelocity stars (HVSs) to set the most precise constraint on the Galactic mass within 80 kpc.

Brown et al. (2010) present a sample of 910 late B-stars and early A-stars in the Galactic halo. Their luminosity, and therefore distance, depends on whether these stars are BHB stars or main-sequence blue stragglers with similar effective temperature and surface gravity. The ambiguous nature of the stars is especially problematic at redder colors, u − g>0.6, where the luminosity can differ by a factor of 5, as demonstrated in Figure 4 in Brown et al. (2010). This bimodal distribution of distance does not have a well-defined average value and therefore requires statistical sampling. Brown et al. (2010) create 100 Monte Carlo realizations of each star to sample the color and metallicity distributions and to derive the distributions of luminosity and distance. We use this full Monte Carlo catalog of distances in our analysis, while retaining the observed values of radial velocity and its uncertainty.

We use the spherical Jeans equation (Binney & Tremaine 2008) to solve for the circular velocity V_c given the radial velocity dispersion σ_r of tracer particles:

where M(r) is the enclosed total mass within radius r, ρ_tr is the density of tracer population, and β is the anisotropy parameter. This equation assumes that the mass distribution M(r) is static and spherically symmetric.

In our case the tracers are halo stars, which follow a steep density profile with negative logarithmic slope (see below)

The velocity dispersion slope and the anisotropy parameter are typically less than unity, and therefore, γ_tr dominates the last factor in Equation (1). The circular velocity profile depends mainly on σ_r and γ_tr.

For ease of comparison with mass profiles expected from cosmological simulations, we fit a power-law relation to the radial velocity dispersion, normalized at radius r₀:

We apply a maximum likelihood (ML) method, similar to that described in the Appendix of Gnedin et al. (2007), to maximize the probability of the model fit given the observations:

We assume that each radial velocity v_r,i is drawn from a Gaussian distribution with zero mean and with combined variance of the intrinsic dispersion at its location, σ_r(r_i), and the measurement uncertainty σ_v,i. The ML method uses all available observational information without binning the data.

The value of the velocity dispersion is sensitive to the presence of outliers in the sample. In order to remove them, we impose a cut |v_r,i| ⩽ V_esc with the escape velocity calculated self-consistently from the derived mass profile. In the "most-likely" case, this process removes only 6 out of 558 stars located at r>25 kpc.

The escape velocity at a given radius r measures the total gravitational potential, a sum of the two terms interior and exterior to r:

We evaluate the second term by adopting a power-law profile for the total matter density ρ(r) ∝ r^−2−α, out to a maximum radius r_out. In the outer Galactic halo we expect 0 < α < 1. Then, the escape velocity is

Ignoring the boundary condition at r_out, we obtain

Watkins et al. (2010) show that most halo tracers reside in the radial range where α ≈ 0.5, and thus V_esc(r) ≈ 2 V_c(r). As the resulting circular velocity profile is close to flat, we adopt this simple relation in our analysis. The ML fit converges to the same value regardless of the initial guess for V_esc. In the inner regions of the Galaxy, lower α leads to higher escape velocity, consistent with the estimate from the RAVE survey of V_esc ≈ 500–600 km s⁻¹ at the solar location (Smith et al. 2007).

The density profile of RR Lyrae halo stars shows a break around 25–30 kpc (Watkins et al. 2009; Sesar et al. 2010), with shallow inner slope (γ_tr ≈ 3) and steep outer slope (γ_tr ≈ 4.5). The stellar halo is also less flattened in the outer regions. In order to obtain a robust measure of the mass profile, we restrict our analysis only to radii r>25 kpc. Figure 1 shows that at these radii the density profiles of both BHB stars and RR Lyrae are consistent with our fiducial value γ_tr ≈ 4.

Zoom In Zoom Out Reset image size

We consider the following parameter ranges: 3.5 ⩽ γ_tr ⩽ 4.5 and 0 ⩽ β ⩽ 0.5. For each parameter set, we use the full Monte Carlo catalog and self-consistently remove outliers using the value of the escape velocity resulting from the velocity dispersion profile. This procedure gives simultaneous estimates for both σ_r(r) and V_c(r).

The radial velocity dispersion profile declines remarkably little with radius. For the power-law fit (Equation (3)), we obtain σ₀ ≈ 114 km s⁻¹ and γ_σ ≈ 0.1 for r₀ = 40 kpc. This choice of r₀ minimizes the error of σ₀. We do not quote formal errors of the fit because the uncertainty is dominated by the unknown parameters γ_tr and β, as we discuss below. Table 1 summarizes the parameters of the best fit and the most extreme allowed models.

Brown et al. (2010) consider a linear fit to the same data and obtain σ_r ≈ 120–0.3 r, where r is in kpc and σ_r is in km s⁻¹. By linearly expanding Equation (3) near r₀, we derive a similar expression in the range 25 kpc < r < 80 kpc: σ_r ≈ 120–0.22 r.

Xue et al. (2008) fit an exponential function to their SDSS DR6 sample and obtain σ_r = 111 exp(−r/354 kpc) km s⁻¹. A linear approximation gives σ_r ≈ 111–0.31 r, which has a similar slope and slightly lower normalization than the HVS sample. Since our dispersion profile is calculated from the Monte Carlo catalog, we similarly bootstrap the Xue et al. (2008) sample 30 times and obtain a distance distribution from which we calculate the dispersion plotted in Figure 2.

Zoom In Zoom Out Reset image size

The profiles from both samples are less steep than that derived by Battaglia et al. (2005): σ_r ≈ 132–0.6 r.

Figure 2 shows the velocity dispersion profiles for all three samples. They appear consistent with each other within the errors, and the differences in the derived best-fit parameters can be attributed to small sample sizes. Outside 25 kpc, the Battaglia, Xue, and Brown samples contain 80, 741, and 558 objects, respectively (using the distance distribution from the Monte Carlo catalogs). Outside 50 kpc, the three samples contain 24, 76, and 163 objects, respectively. Note that our sample more than doubles the number of distant stars and thus presents the most accurate current measurement of σ_r. Individual velocity errors and resampling of the distance distribution affect the value of σ_r by less than 3%. The streaming velocity 〈v_r〉 is small everywhere except at the outermost radii. Overall, the data indicate that the velocity dispersion varies little with radius out to 80 kpc.

We can now combine Equations (1) and (3) to calculate the circular velocity

This estimate is degenerate with respect to the density and anisotropy of the tracer population. Lower γ_tr and higher β result in the lower mass estimate and vice versa. At the same time, steeper density slope can balance stronger anisotropy. Based on our derivation of the tracer density in Figure 1, we believe that the slope is constrained to be between 3.5 and 4.5, with the most likely value of γ_tr ≈ 4. The anisotropy parameter is not directly known, but we can take the predictions of cosmological simulations of galaxy formation as a guide. The centers of halos are close to being isotropic. In the outer parts, the orbits of dark matter particles and satellite halos become more radially biased with distance, reaching β ≈ 0.5 outside the peak of the circular velocity curve (e.g., Diemand et al. 2007; Navarro et al. 2010). Thus, a full range of values 0 ⩽ β ⩽ 0.5 is possible, while the most likely value in our radial range is β ≈ 0.4.

Figure 2 illustrates the allowed range of circular velocity. We take the best fit (γ_tr = 4, β = 0.4) and consider the full range of variation of either parameter while keeping the other fixed. The overlap of the two shaded regions gives the most probable values of V_c. The uncertainty due to either parameter is systematic in nature and comparable in magnitude. The value of the circular velocity at 80 kpc is uncertain at least by 10% and the value of the enclosed mass at least by 20%. The most likely values are V_c(80) = 193 km s⁻¹, M(80) = 6.9 × 10¹¹ M_☉.

For comparison, we also use an independent robust mass estimator proposed by Watkins et al. (2010). The mass within radius r is given by the following average of radial velocities of all objects inside r:

where (−2 − α) is again the logarithmic slope of the density profile in the range of radii probed by the data. This estimator also depends on γ_tr and β, similarly to Equation (8). In fact, if these two parameters and α are constant with radius, then α = 2γ_σ. The estimator is different in that it takes a single average of all radial velocities instead of fitting the dispersion profile.

Three vertical lines in Figure 2 show the applications of the robust estimator at 40, 60, and 80 kpc. They indicate a slightly steeper decline of V_c with radius, but still fully consistent with the results derived from the dispersion. The most likely value at 80 kpc is $V_{c,\rm est}\, (80) \approx 190\ \mathrm{km\ s}^{-1}$.

Figure 3 compares our derived circular velocity with other independent observational constraints. Dotted lines show the extrapolation of V_c outside the range of our data and should be treated with caution. Watkins et al. (2010) estimate the mass of the Galaxy at r = 100 kpc using the sample of all known satellite galaxies. Depending on the anisotropy and the inclusion of particular galaxies, the mass can vary between 3.3 and 13.8 × 10¹¹ M_☉ corresponding to the circular velocity range from 119 to 244 km s⁻¹. Our results favor the upper half of this range.

Zoom In Zoom Out Reset image size

Within 60 kpc, Xue et al. (2008) estimate V_c = 170 ± 15 km s⁻¹ by matching observed velocity dispersion to the motion of particles in two hydrodynamic simulations of galaxy formation. Within 50 kpc, Wilkinson & Evans (1999) obtain M(50) = 5.4^+0.2_−3.6 × 10¹¹ M_☉ by modeling the distribution function of observed velocities. For clarity, we do not show similar results obtained by Kochanek (1996), M(50) = 4.9^+1.1_−1.1 × 10¹¹ M_☉, and Sakamoto et al. (2003), M(50) = 5.5^+0.1_−0.4 × 10¹¹ M_☉.

The circular velocity at the solar circle, V₀, is constrained better but it scales with the distance to the Galactic center, R₀. McMillan & Binney (2010) analyze the motion of masers in star-forming regions throughout the Galaxy and find that the best-constrained parameter is the ratio V₀/R₀ = 30  ±   2 km s⁻¹ kpc⁻¹. The distance R₀ can vary from 6.7 to 8.9 kpc, depending on model assumptions. We plot this constraint on V₀ as a diagonal line in Figure 3. This line matches the extrapolation of our fits perfectly. However, this result should not be considered as an improvement on the value of V₀, because we do not expect V_c(r) to remain a power law with fixed slope at such small radii.

For illustration, we also show an example of a three-component model of the Galactic potential, similar to the cosmologically motivated model by Klypin et al. (2002). We assume an axisymmetric exponential disk with the mass 5 × 10¹⁰ M_☉ and scale length 3 kpc, and a compact bulge with the mass 5 × 109 M_☉ (Binney & Tremaine 2008). We include a dark matter halo represented by a Navarro–Frenk–White (NFW) profile with the scale radius r_s = 20 kpc and the mass that, along with the disk and the bulge, gives the total virial mass M_vir = 1.6 × 10¹² M_☉. This model lies in the middle of our inferred interval of circular velocity and can be taken as a good first approximation to the Galactic mass distribution. Note however that the uncertainty in the value of the virial mass is at least 20% and possibly larger.

In Figure 4, we further illustrate the dependence of the derived value of the circular velocity at 80 kpc on the parameters γ_tr and β. For each combination of the parameters, we run an ML fit to the velocity dispersion profile with self-consistent removal of unbound stars. The plot shows, however, that contours of V_c(80) lie close to the lines γ_tr − 2β = const, as expected from Equation (8) for a constant γ_σ. Larger β gives smaller velocity, but according to cosmological simulations it should not exceed 0.5 and be near 0.4. However, the estimate of V_c would be higher if the tracers had tangential anisotropy, β < 0. A density slope γ_tr ≈ 4 is most likely, as it fits both the BHB and RR Lyrae samples.

Zoom In Zoom Out Reset image size

Using ML analysis of a new sample of radial velocities of distant halo stars, we infer that their radial velocity dispersion profile declines little with distance from the Galactic center in the range 25 kpc < r < 80 kpc: σ(r) = 111(r/40 kpc)^−0.08 km s⁻¹. Spherical Jeans modeling indicates that the circular velocity profile V_c(r) also falls only slightly over the same radial range and reaches between 175 and 231 km s⁻¹ at 80 kpc. The corresponding enclosed mass M(80) is between 5.7 × 10¹¹ M_☉ and 1.0 × 10¹² M_☉. A three-component model for the baryon and dark matter mass distribution gives the total virial mass of the Galaxy M_vir = (1.6 ± 0.3) × 10¹² M_☉ at the virial radius R_vir = 300 kpc.

Our inferred mass of the Galaxy is higher than that obtained by Battaglia et al. (2005; M_vir ≈ 0.8 × 10¹² M_☉) and Xue et al. (2008; M_vir ≈ 1.0 × 10¹² M_☉) based on the modeling of their radial velocity data sets. Our HVS sample contains more objects at r>40 kpc than the Battaglia et al. and Xue et al. data sets. Thus, we have a stronger constraint on the shallow slope of the velocity dispersion profile and we derive a correspondingly larger mass. Our inferred mass is consistent with the larger scale measurement by Li & White (2008) based on the Andromeda Milky Way timing argument, M_vir ≈ 2.4 × 10¹² M_☉. The implied dynamical mass-to-light ratio of the Galaxy, M_vir/L_V ≈ 50 in solar units, is also consistent with galaxy–galaxy weak-lensing measurements by Mandelbaum et al. (2006), galaxy kinematics modeling by More et al. (2010), and halo abundance matching modeling by Moster et al. (2010).

Our sample of radial velocities is large enough that the biggest uncertainty in the mass estimate is not statistical but systematic. Within the framework of spherical Jeans modeling, the uncertainty is dominated by the density slope and anisotropy of the tracer population. These parameters could be better constrained by future all-sky surveys of halo BHB stars. Deeper surveys that target more distant stars at r ≳ 100 kpc would be similarly dominated by uncertainty over the underlying distribution of the tracers.

The validity of spherical Jeans modeling is also limited by the presence of structure in the distribution of halo stars. Galactic stellar halo contains remnants of disrupted satellite galaxies, some of which are still detectable as tidal streams. Stars at ∼100 kpc from the Galactic center may not have had enough dynamical times to reach dynamical equilibrium, further limiting the application of equilibrium modeling. A first step in the direction of circumventing these systematics was taken by Xue et al. (2008), who modeled the motion of particles in realistic simulated halos. Extension of such analysis to many different halo realizations using large samples of observed velocities may reduce the uncertainty over the global mass distribution in the Galaxy.

We thank Hans-Walter Rix for clarifying discussions and the organizers of the Sixth Harvard-Smithsonian Sackler Conference, where this work was completed. O.G. is supported in part by the NSF Grant AST-0708087.