VELOCITY DISPERSION PROFILE OF THE MILKY WAY HALO

All photometry comes from SDSS Data Release 6 (Adelman-McCarthy et al. 2008). We use uber-calibrated point-spread function (PSF) magnitudes, and correct the magnitudes and colors for reddening following Schlegel et al. (1998).

The HVS survey target selection emphasizes outliers in the halo population: we target stars redder in (u − g)₀ than known white dwarfs and bluer in (g − r)₀ than known BHB stars (Brown et al. 2006b). This color cut through the stellar population, illustrated in Figure 1, allows us to detect HVSs efficiently. The majority of targets, however, are normal halo stars.

Zoom In Zoom Out Reset image size

Our target selection includes stars with 17 < g₀ < 19.5 in the range −0.39 < (g − r)₀ < −0.25 (Brown et al. 2007c) and fainter stars with 19 < g₀ < 20.5 over a broader color range −0.40 < (g − r)₀ < −0.20 (Brown et al. 2009a). We also include eight confirmed BHB stars from Brown et al. (2005a) with 19.5 < g₀ < 20.25 and −0.3 < (g − r)₀ < −0.1. All 910 targets are located in the SDSS DR6 footprint and have an average surface density on the sky of 0.12 deg⁻².

We obtained spectroscopic observations at the 6.5 m MMT telescope with the Blue Channel Spectrograph. We operated the spectrograph with the 832 line mm⁻¹ grating in second order, providing wavelength coverage 3650–4500 Å and a spectral resolution of 1.2 Å. We obtained all observations at the parallactic angle, with a comparison lamp exposure for every survey object.

We processed the data using IRAF⁴ in the standard way. We measure radial velocities by cross-correlating the observations with radial velocity standards (Fekel 1999) using the package RVSAO (Kurtz & Mink 1998). The average radial velocity uncertainty of the stars is ±12 km s⁻¹.

All velocities discussed here are in the Galactocentric rest frame, indicated v_rf. Given the outer halo location of the stars, we note that the observed radial velocities are almost purely (>85%) radial in the Galactocentric frame. We transform heliocentric velocities (v_helio) into Galactocentric rest-frame velocities assuming a circular velocity of 220 km s⁻¹ and a solar motion of (U, V, W) = (10, 5.2, 7.2) km s⁻¹ (Dehnen & Binney 1998):

$$\begin{eqnarray} v_{\rm rf} &=& v_{\rm helio} + 220\sin {l}\cos {b} + (10\cos {l}\cos {b}\nonumber\\ &&+ 5.2\sin {l}\cos {b} + 7.2\sin {b}). \end{eqnarray} \tag{ 1 }$$

Reid et al. (2009) argue for a larger circular velocity of 250 km s⁻¹ based on trigonometric parallaxes to star formation regions in the disk. We test using a circular velocity of 250 km s⁻¹ and find statistically identical velocity dispersion profiles. This insensitivity to circular velocity arises because the stars in our high latitude survey have a mean value of |sin lcos b| = 0.32; changing the Sun's circular velocity by 30 km s⁻¹ results in a ±10 km s⁻¹ change in the stars' rest-frame velocities. 10 km s⁻¹ is smaller than our measurement error and an order of magnitude smaller than the velocity dispersion of the stars. McMillan & Binney (2009) show that the Reid et al. (2009) data are consistent, at the 1σ level, with the canonical circular velocity of 220 km s⁻¹. Thus, we use 220 km s⁻¹ here.

Figure 2 plots the resulting Galactic rest-frame velocity distribution of the 910 halo stars. For reference, we also draw a Gaussian with zero mean and 106 km s⁻¹ dispersion (dotted line). The observations reveal a significant asymmetry of outliers in the tails of the velocity distribution. There are no stars with v_rf < −300 km s⁻¹ but 18 stars with v_rf > + 300 km s⁻¹(the HVSs). We address this issue in Section 3.

Zoom In Zoom Out Reset image size

The strongest metal line in our spectra is the 3933 Å Ca ii K line. Unfortunately, at the effective temperatures sampled by the survey, 10,000 ⩽ T_eff ⩽ 15,000 K, the equivalent width of Ca ii K is small (Wilhelm et al. 1999a). Thus, Ca ii K provides poor leverage on the metallicity of our stars. Metallicity is better determined for redder (cooler) BHB stars, for example, the BHB survey of Brown et al. (2008a) and the BHB sample from the SDSS spectroscopic survey (Xue et al. 2008). The metallicity distributions of these two BHB samples are plotted in Figure 3 and are similar in shape. However, the outer halo sample of Xue et al. (2008) is about 0.2 dex more metal poor than the inner halo sample of Brown et al. (2008a). The mean metallicity of the Xue et al. (2008) BHB stars with g₀ > 17 is [Fe/H]_Ca = −1.9.

Zoom In Zoom Out Reset image size

To make luminosity estimates, we assume that our survey stars have the metallicity distribution function of Xue et al. (2008). This assumption is reasonable given the very similar sky coverage and penetration into the halo of the two surveys: the stars occupy similar regions of the Milky Way halo. The reddest stars in our sample, where we can estimate metallicity, are metal poor (Brown et al. 2006b), consistent with the metallicity distribution function of Xue et al. (2008).

Although the stars in our survey have the spectral types of late-B- and early-A-type stars, their nature is ambiguous. The old stellar population of the halo contains both evolved BHB stars and main-sequence blue stragglers. Halo surveys consistently find that ∼50% of A-type stars in the field are blue stragglers (Norris & Hawkins 1991; Kinman et al. 1994; Preston et al. 1994; Wilhelm et al. 1999b; Clewley et al. 2002, 2004; Brown et al. 2003, 2005b, 2008a; Xue et al. 2008). This result is problematic because BHB stars and blue stragglers can have very different luminosities (Figure 1). Spectroscopic measures of surface gravity can discriminate the evolutionary state of the stars (Kinman et al. 1994; Wilhelm et al. 1999a; Clewley et al. 2002, 2004). Unfortunately, surface gravity measures fail for our sample because BHB and main-sequence stars have nearly identical surface gravities at the effective temperatures of our stars.

We target stars so blue, however, that we largely exclude the possibility of blue stragglers. A ≃0.75 M_☉ star with [Fe/H] = −1.9 has a main-sequence lifetime of a Hubble time (Girardi et al. 2004). A field blue straggler cannot plausibly have more than twice the mass of a main-sequence turnoff star. In our sample, 36% of the stars are bluer than both (g − r)₀ = −0.27 and (u − g)₀ = 0.87, the color of a 1.5 M_☉ star with [Fe/H] = −1.9 (Girardi et al. 2004). Thus, the nature of many of the stars is clear: they are hot BHB stars.

Stars bluer than a 1.5 M_☉ star are BHB; stars redder than a 1.5 M_☉ star are equally likely to be BHB stars or blue stragglers. We base this conclusion on Xue et al. (2008), who observe a BHB fraction of 47% at the red end of our sample. Overall, our sample is 74% BHB stars and 26% blue stragglers.

We estimate luminosity from stellar evolution tracks. Girardi et al. (2002, 2004) provide main-sequence tracks in the SDSS passbands for metallicities ranging from solar to [Fe/H] = −2.3, but do not include the horizontal branch. Dotter et al. (2007, 2008) provide BHB tracks in the SDSS passbands over the same range of metallicities. To simplify the luminosity estimates, we fit low-order polynomials to the tracks as a function of both color and metallicity.

We estimate luminosity two ways, using a star's (u − g)₀ and (g − r)₀ color. Our final luminosity is the weighted average of the two estimates. We weight by the slope of the observed (u − g)₀ versus (g − r)₀ distribution (Figure 1). This weighting is important because luminosity is most sensitive to (u − g)₀ at the blue end of our sample and to (g − r)₀ at the red end of our sample. The weighting favors luminosities derived from (u − g)₀ for (g − r)₀ < −0.24 and luminosities derived from (g − r)₀ for (g − r)₀ > − 0.24. Figure 1 shows the resulting lines of constant luminosity for a star with [Fe/H] = −1.9 from the Girardi et al. (2004) main-sequence tracks (in blue) and from the Dotter et al. (2008) BHB tracks (in red).

The accuracy of the luminosity estimate depends on the accuracy of the stellar evolution models. For example, alpha-enhanced tracks (appropriate for the old stellar halo populations) are systematically bluer than solar-scaled tracks for both main-sequence and BHB stars (Lee et al. 2009).

The precision of the luminosity estimate is more robust because we apply the same tracks to all stars. The precision depends on observational uncertainties and is straightforward to quantify. A g₀ = 19 star at the median depth of our survey has uncertainties σ_{(u − g)} = 0.05 and σ_{(g − r)} = 0.03 in color and σ_v = 12 km s⁻¹ in velocity. We assume that metallicity is randomly drawn from the Xue et al. (2008) metallicity distribution function (Figure 3). We then propagate these uncertainties through the main-sequence and BHB tracks to visualize the resulting distribution of luminosity estimates.

Figure 4 plots the distribution of luminosity and velocity (centered at zero) for two fiducial stars in the bluest quartile (left panel) and the reddest quartile (right panel) of the survey. Stars in the bluest quartile have luminosities precise to 10% because main-sequence and BHB stars have essentially identical luminosities at these colors (Figure 1). Stars in the reddest quartile, however, have bimodal luminosity distributions with a factor of 4 spread in luminosity (a factor of 2 in distance). This degeneracy is broken for confirmed BHB stars: 33 of our stars are listed as BHB stars in Xue et al. (2008) and 8 are BHB stars from our earliest sample (Brown et al. 2005a). For these confirmed BHB stars, we use only the BHB luminosity to calculate distance.

Zoom In Zoom Out Reset image size

Three of the stars with v_rf > + 400 km s⁻¹ in this survey are confirmed young main-sequence B stars (Fuentes et al. 2006; Bonanos et al. 2008; López-Morales & Bonanos 2008; Przybilla et al. 2008a, 2008b). For the purposes of this paper, however, we treat these HVSs the same as the other stars in our sample: we estimate distances to the HVSs as if they were halo BHB stars or metal-poor blue stragglers. We clip the HVSs from the sample when calculating the velocity dispersion in Section 3.

We calculate distances from the observed apparent magnitudes and estimated luminosities:

$$\begin{equation} d = 10^{(g_0 - M_g)/5-2}\;{\rm kpc}, \end{equation} \tag{ 2 }$$

where d is the heliocentric distance in kpc, g₀ is the apparent magnitude corrected for extinction, and M_g is the absolute magnitude estimated above. We convert heliocentric distance d to Galactocentric distance R with the Sun at R = 8 kpc.

Figure 5 plots the resulting distance distribution of the survey stars compared with the Xue et al. (2008) sample of BHB stars. The Xue et al. (2008) sample contains more than twice as many stars as our sample; however, our sample is deeper and contains twice as many stars with R ≳ 50 kpc. Neither sample fairly measures the density profile of the halo. The Xue et al. (2008) sample is incomplete in color, magnitude depth, and spatial coverage. Our sample is complete in all dimensions, but we find an artificially shallow density profile—consistent with N(R) ∝ R⁻² over the range 20 < R < 70 kpc (Figure 5)—because we observe over a broader color range starting at g₀ > 19 mag. A shallow density profile works to our advantage for the velocity dispersion analysis, however, because our stars sample greater distances at a greater relative density.

Zoom In Zoom Out Reset image size

We use a Monte Carlo approach to model the distance to each star. We assume that photometric and velocity errors are Gaussian distributed, and that the underlying metallicity distribution is that of Xue et al. (2008). We then derive the luminosity for each star by randomly drawing its color and metallicity from the observed distributions and comparing them to the main-sequence and BHB tracks (Figure 1) described above. We sample the distributions 100 times per star. Using a Monte Carlo approach allows us to account for the non-Gaussian distribution of luminosity estimates unique to each star.

The resulting Monte Carlo catalog drawn from the observations produces the "cloud" of velocities and distances shown in Figure 6; the distribution of points reflects the uncertainties in the velocity measurements and the distance estimates of the stars.

Zoom In Zoom Out Reset image size

We calculate the velocity dispersion profile by grouping stars into five or six unique bins in distance. The bins have sizes of 0.33R, chosen so that the bins contain at least 70 stars. We require 50 stars to obtain a dispersion with a formal statistical uncertainty of 10%. Bins with smaller occupation have dispersions systematically biased low resulting from small number statistics. Given the distance distribution of our sample, the occupation requirement constrains our velocity dispersion measurements to the region 15 < R < 75 kpc. Stars with well determined luminosities contribute their full weight to a single distance bin; stars with poorly constrained luminosities contribute less weight distributed over multiple bins.

We use bootstrap re-sampling to calculate the uncertainty in the velocity dispersion measurement. Our procedure is to draw random sets of 910 stars, with replacement, from the Monte Carlo catalog and re-calculate the velocity dispersion. We re-sample 10,000 times; the uncertainty is the standard deviation of the velocity dispersions.

Our survey was designed to find unbound HVSs. Thus, the first step in calculating the velocity dispersion profile is to clip the unbound stars from the sample. Defining an unbound star is difficult, however. The best estimate of the Galactic escape velocity is 550 ± 50 km s⁻¹ at the solar circle, based on the highest velocity stars in the solar neighborhood (Smith et al. 2007). The escape velocity of the outer halo is more poorly constrained; extrapolating v_esc(R) to large R requires a potential model.

Here, we use the Kenyon et al. (2008) potential model that is tied to observed Milky Way mass measurements. To establish v_esc(R) in this model, we drop a test particle at rest from R = 500 kpc (i.e., half-way to M31) and calculate its trajectory down to R = 0 kpc. The Kenyon et al. (2008) model predicts v_esc = 585 km s⁻¹ at R = 8 kpc in this definition of escape velocity, in good agreement with the Smith et al. (2007) observation. We fit a low-order polynomial to the calculated velocity as a function of distance and find

$$\begin{equation} v_{\rm esc}(R) = -2.30 \, \times \, 10^{-4} R^3 \, + \, 0.0588 R^2 \, - \, 6.62 R \, + \, Z \ {\rm km\,s^{-1}}, \end{equation} \tag{ 3 }$$

where Z = 619 km s⁻¹, valid for 15 < R < 100 kpc.

To avoid imposing an arbitrary mass on our velocity dispersion measurement, however, we re-normalize Equation (3) to the observed envelope of negative velocity stars in our sample. Figure 6 shows the v_esc(R) relation that we use, with Z = 500 km s⁻¹ (Equation (3)). Our assumption in changing Z to 500 km s⁻¹ is that the velocity distribution of halo stars extends up to the escape velocity, and that the most negative velocity stars—the stars falling in from the largest distances—provide the most robust measure of escape velocity in our data set. Perets et al. (2009) discuss this further in the context of HVSs. We test the effects of our choice of v_esc(R) by using a non-parametric approach to calculate velocity dispersion in the next section. In this section, we calculate the velocity dispersion profile by using the shape of the v_esc(R) relation to provide a physically motivated means of clipping velocity outliers.

We plot the resulting velocity dispersion of stars with velocities less than v_esc(R) for Z = 500 km s⁻¹ in the lower panel of Figure 6. Our Monte Carlo approach allows outliers that might otherwise be clipped to contribute a weight appropriate to their measurement uncertainties. A linear least-squares fit to the velocity dispersion profile finds a declining velocity dispersion,

$$\begin{equation} \sigma _v = (-0.30 \pm 0.10)R + (120\pm 4.6)\;{\rm km\;s^{-1}}, \end{equation} \tag{ 4 }$$

valid over 15 < R < 75 kpc. A higher-order fit does not significantly improve the residuals. The linear fit has a standard deviation of 4.8 km s⁻¹ and a reduced χ² of 0.7. A fixed velocity dispersion, by comparison, has a standard deviation of 7.6 km s⁻¹ and a reduced χ² of 1.5, poorer than the linear fit but not statistically inconsistent.

A non-parametric approach to determining the velocity dispersion profile leads to robust results based on fewer a priori assumptions. Here, we use the caustic technique. Because this technique is new to the stellar halo literature, we begin with an introduction before discussing the mechanics and results. Diaferio (2009) provides a recent review of the technique.

The caustic technique was originally developed to measure galaxy cluster mass profiles (Diaferio & Geller 1997; Diaferio 1999). Caustics are essentially escape velocity curves. The caustic technique is based on the distribution of objects in phase space: the line-of-sight velocity v versus projected distance r from the cluster center. In this plane, cluster members distribute in a characteristic trumpet shape with upper and lower borders indicating the escape velocity from the system. The amplitude ${\cal A}(r)$ of this distribution, which decreases with r, is a direct measure of the escape velocity from the system, independent of its dynamical state. The caustic technique only assumes that the spatial distribution of tracers is spherically symmetric.

The caustic technique can be applied to any self-gravitating system. The technique is non-parametric and does not require that the tracers of the potential be in virial equilibrium. It is an effective tool for identifying bound members of the system defined by the upper and lower caustics (sharp declines in phase space density of tracers). For example, Serra et al. (2009) apply the caustic technique to five dwarf spheroidals of the Milky Way.

Applying the caustic technique to the halo star sample requires a slight modification of the procedure used for galaxy clusters. Traditionally, the technique arranges sample objects in a binary tree according to their pairwise projected binding energies (Diaferio 1999). By walking along the main branch of the tree, the technique determines the velocity dispersion which locates the caustic in the velocity diagram. This step does not apply to the Milky Way halo; the stars are not obviously clustered on the sky and the estimate of the pairwise binding energy is inappropriate. Thus, we apply the technique directly to the phase space diagram shown in Figure 6.

The caustics are the curves satisfying the equation f_q(r, v) = κ, where f_q(r, v) is the distribution in the phase space diagram, and κ is the root of the equation 〈v²_esc〉_κ,R = 4σ². The function

$$\begin{equation} \big\langle v_{\rm esc}^2\big\rangle _{\kappa,R}=\int _0^R{\cal A}_\kappa ^2(r)\varphi (r)dr/ \int _0^R\varphi (r)dr \end{equation} \tag{ 5 }$$

is the mean caustic amplitude within R, φ(r) = ∫f_q(r, v)dv. R is not a free parameter in the standard application of the technique, but here we chose R = 30 kpc. Our results are totally insensitive to this parameter: varying R in the range 20–100 kpc varies the final number of the halo members by at most five and does not change the velocity dispersion profile.

We use an iterative approach to estimate the velocity dispersion profile. First, we compute the velocity dispersion of the total sample and locate the caustics that enable some interloper removal. We then compute a new velocity dispersion with the member stars, locate new caustics, and remove further interlopers. We proceed until no star is identified as an interloper. This procedure requires only four steps and on average yields 838 (of 910) final halo members and an overall velocity dispersion of 99 km s⁻¹. The final set of caustics is drawn with the solid line in the upper panel of Figure 6.

The open triangles in the lower panel of Figure 6 show the caustic velocity dispersion profile. It is visually apparent that the caustic technique measures the velocity dispersion from the well-sampled core of the velocity distribution. This approach results in a reliable velocity dispersion profile, but yields a velocity dispersion that is systematically 10% smaller than found by our parametric v_esc(R) method. We conclude that clipping velocity outliers with the caustic technique changes the observed amplitude, but not the observed slope, of the velocity dispersion profile.

A linear least-squares fit to the caustic velocity dispersion profile finds a declining velocity dispersion,

$$\begin{equation} \sigma _v=(-0.45 \pm 0.16) R + (110\pm 6.0)\;{\rm km\;s^{-1}}, \end{equation} \tag{ 6 }$$

valid over 16 < R < 64 kpc. The linear fit has a standard deviation of 4.8 km s⁻¹ and a reduced χ² of 1.0. A fixed velocity dispersion, by comparison, has a standard deviation of 7.9 km s⁻¹ and a reduced χ² of 2.8. In this case, a constant velocity dispersion is a significantly poorer fit to the data.

There are at least two contaminants that may systematically affect our velocity dispersion profile: binary stars, which increase the velocity dispersion, and disk stars (i.e., white dwarfs), which decrease the velocity dispersion. We investigate how these contaminants may alter the observed velocity dispersion profile.

Binary stars are unlikely to change the velocity dispersion profile. Both BHB stars and 1.5 M_☉ blue stragglers have stellar radii around 2.5 R_☉. Equal mass pairs of such stars must have semimajor axes of at least 6.5 R_☉ to avoid Roche lobe overflow. Thus, the most compact possible binary system has velocity semiamplitude of 100 km s⁻¹. Assuming that half of the targets are binaries with a lognormal distribution of semimajor axes, we expect 30 binaries in our data set with vsin i > 50 km s⁻¹. This estimate is generous given that BHB stars have recently evolved through the red giant phase and are thus unlikely to have close companions. In any case, binaries will be observed at a random orbital phase and binned with 70–200 other stars to compute a velocity dispersion. We propagate our simulated distribution of vsin i's through our Monte Carlo catalog and find a negligible change in the velocity dispersion profile; the uncertainty in the velocity dispersion slope remains 30%.

Disk stars pose a greater threat to the velocity dispersion profile. Disk stars have a systematically lower velocity dispersion than halo stars and may also have a non-zero mean velocity because of the longitude dependence of the Sun's circular velocity correction. The SDSS imaging survey from which we draw our candidates does not uniformly survey the sky across all longitudes.

Our greatest concern is white dwarfs, which are intrinsically faint and overwhelmingly appear at faint magnitudes in our survey. Inserting white dwarfs into our halo sample thus reduces the velocity dispersion observed at large distances. Observationally, 15% of our survey targets are white dwarfs, mostly found at low latitudes. If we insert a 1% white dwarf contamination into our Monte Carlo catalog, the resulting velocity dispersion profile steepens by 10%. Fortunately, white dwarfs are readily identified by their broad Balmer lines (see Kilic et al. 2007a). Visual inspection of our spectra reveals no white dwarf contaminants.

Binaries also have observational constraints. We obtained repeat observations for every star with |v_rf|>300 km s⁻¹ and found only one star, a star near −300 km s⁻¹, to exhibit radial velocity variation. We thus excluded this star from the halo sample. We conclude that binaries and white dwarfs are unlikely to significantly influence the velocity dispersion profile measured from our halo star sample.

There are few measurements of the Milky Way velocity dispersion profile at large R. The most comparable measurements are Battaglia et al. (2005) and Xue et al. (2008). Battaglia et al. (2005) measure the velocity dispersion profile based on a sample of 9 satellite galaxies, 44 globular clusters, 58 red giants, and 130 BHB stars that span 10 < R < 140 kpc. Xue et al. (2008) measure the velocity dispersion profile based on a sample of 2401 BHB stars that span 5 < R < 60 kpc. With the exception of a handful of shared BHB stars, the samples are independent of one another.

Figure 7 compares the velocity dispersion profile and radial number distribution of our sample with those of Battaglia et al. (2005) and Xue et al. (2008). Beyond R > 50 kpc, our sample contains a factor of 8 more stars than Battaglia et al. (2005) and a factor of 2 more stars than Xue et al. (2008). Remarkably, the velocity dispersions measured by the three samples are consistent at the 1.5σ level. The three samples also observe a statistically similar decline in velocity dispersion with distance. The major difference between the samples is the significance of the observed decline in velocity dispersion with distance.

Zoom In Zoom Out Reset image size

A linear least-squares fit to the Xue et al. (2008) data gives a −0.15 ± 0.14 km s⁻¹ kpc⁻¹ decline in velocity dispersion over 10 < R < 60 kpc. The linear fit has a standard deviation of 6.3 km s⁻¹ and a reduced χ² of 1.5. A fixed velocity dispersion has a very similar standard deviation of 6.4 km s⁻¹ and reduced χ² of 1.9. Thus, the Xue et al. (2008) measurements cannot formally discriminate between a constant velocity dispersion and a declining velocity dispersion.

The Battaglia et al. (2005) data, on the other hand, exhibit a steeper −0.58 ±  0.11 km s⁻¹ kpc⁻¹ decline in velocity dispersion over 10 < R < 100 kpc. We exclude their last bin because it contains only three objects. The linear fit has a standard deviation of 8.6 km s⁻¹ and a reduced χ² of 0.85.

The weighted average of all the samples yields a −0.38 ± 0.12 km s⁻¹ kpc⁻¹ decline in velocity dispersion over 10 < R < 100 kpc. We obtain the same result if we average only our own two velocity dispersion measurements. The Milky Way velocity dispersion profile is not linear, of course, but three independent sets of observations are in statistical agreement: the Milky Way radial velocity dispersion drops from σ ≃ 110 km s⁻¹ at R = 15 kpc to σ ≃ 85 km s⁻¹ at R = 80 kpc.