ABSTRACT
We characterize the radial density, metallicity, and flattening profile of the Milky Way's stellar halo, based on the large sample of spectroscopically confirmed giant stars from SDSS/SEGUE-2, spanning galactocentric radii 10 kpc ≤ rGC ≤ 80 kpc. After excising stars that were algorithmically attributed to apparent halo substructure (including the Sagittarius stream), the sample has 1757 K giants, with a typical metallicity precision of 0.2 dex and a mean distance accuracy of 16%. Compared to blue horizontal branch stars or RR Lyrae variables, giants are more readily understood tracers of the overall halo star population, with less bias in age or metallicity. The well-characterized selection function of the sample enables forward modeling of those data, based on ellipsoidal stellar density models, ν*(R, z), with Einasto profiles and (broken) power laws for their radial dependence, combined with a model for the metallicity gradient and the flattening profile. Among models with constant flattening, these data are reasonably well fit by an Einasto profile of n = 3.1 ± 0.5 with an effective radius and a flattening of q = 0.7 ± 0.02, or comparably well by an equally flattened broken power law, with radial slopes of αin = 2.1 ± 0.3 and αout = 3.8 ± 0.1, with a break radius of rbreak = 18 ± 1 kpc; this is largely consistent with earlier work. We find a modest but significant metallicity gradient within the "outer" stellar halo, [Fe/H] decreasing outward. If we allow for a variable flattening , we find the distribution of halo giants to be considerably more flattened at small radii, q(10 kpc) = 0.55 ± 0.02, compared to q(>30 kpc) = 0.8 ± 0.03. Remarkably, the data are then very well fit by a single power law with index of 4.2 ± 0.1 on the variable . In this simple and better-fitting model, there is a break in flattening at ∼20 kpc, instead of a break in the radial density function. While different parameterizations of the radial profile vary in their parameters, their implied density gradient, , is stable along a direction intermediate between major and minor axis; this gradient is crucial in any dynamical modeling that uses halo stars as tracers.
Export citation and abstract BibTeX RIS
1. INTRODUCTION
The Milky Way's extended stellar halo contains only a small fraction (≲1%) of the Galaxy's stars but is an important diagnostic of its formation and dark matter distribution. The position-kinematics-abundance substructure in the stellar halo reflects the Galaxy's formation history, whether halo stars were born in situ or are disrupted satellite debris. By now, individual stars are by far the largest sample of kinematic tracers, with (as opposed to globular clusters or satellite galaxies), and are hence the best tracers to determine the mass profile of Milky Way's dark matter halo in this radial range. It is obvious that good kinematic tracer samples should be large and cover a wide radial range with accurate individual distances. However, beyond this, the spatial distribution of the tracers—in particular their radial profile—must be well understood to use such tracers for dynamical inferences. This point is perhaps clearest when considering the Jeans (1915) equation, even in its simplest, spherical and isotropic version: the tracer density profile, ν*(r), in particular its logarithmic radial derivative, , affects the inferred enclosed mass, M(<r), almost linearly. If we do not know the local power-law exponent to better than, say, 25%, we cannot infer the mass to better than 25% irrespective of the size and quality of the kinematic sample. To a somewhat lesser extent, the inferred mass also depends on the flattening of the tracer population. Yet, at present there is little consensus on the shape and the radial profile of the stellar halo.
The most straightforward way to quantify the stellar halo distribution is via star counts. Because this method requires large samples of well-understood completeness, it is often applied to photometric catalogs. Early studies adopted star counts to analyze globular clusters (Harris 1976), RR Lyrae variables (Hawkins 1984; Wetterer & McGraw 1996; Vivas & Zinn 2006), blue horizontal branch (BHB) stars (Sommer-Larsen 1987), a combination of BHBs and RR Lyraes (Preston et al. 1991), a star sample near the north galactic pole (Soubiran 1993), or K dwarfs (Gould et al. 1998) and found that the stellar halo is well fitted by a single power law (SPL; ) with index α = 3–3.5 and flattening of q = 0.5–1. However, Saha (1985) found that RR Lyrae are well described by a broken power law (BPL) with α ∼ 3 out to 25 kpc, and α ∼ 5 beyond 25 kpc.
These earlier studies were based on a few hundred objects at most. In recent years, with the development of large sky surveys, the sizes of the photometric samples have expanded by a factor of more than 10. Robin et al. (2000) used a wide set of deep star counts in a pencil-beam survey at high and intermediate Galactic latitudes to model the density profile and found the best-fit density profile with a flattening of 0.76 and a power index of 2.44. Siegel et al. (2002) found that 70,000 stars in seven Kapteyn selected areas are consistent with a power-law density with index of 2.75 and flattening of 0.6. However, both Robin et al. (2000) and Siegel et al. (2002) used galaxy models to work out the contamination of star counts by the disk population. This technique is not ideal for the lowest-density stellar component in the Galaxy because any errors in the disk or thick disk will overwhelm the halo results. Morrison et al. (2000) used halo stars near the main-sequence (MS) turnoff from the "SPAGHETTI" survey to map the Galactic halo and found that the halo density law over galactocentric radii of 5–20 kpc and z-heights of 2–15 kpc followed a flattened power-law halo with the flattening of 0.6 and a power index of 3. Bell et al. (2008) used ∼4 million color-selected MS turnoff stars from the fifth data release of the Sloan Digital Sky Survey (SDSS) up to 40 kpc to find a "best-fit" oblateness of the stellar halo of 0.5–0.8, and the density profile of the stellar halo is approximately described by a power law with index of 2–4. Subsequently, Sesar et al. (2011) used 27,544 near-turnoff MS stars out to ∼35 kpc selected from the Canada–France–Hawaii Telescope Legacy Survey to find a flattening of the stellar halo of 0.7 and the density distribution to be consistent with a BPL with an inner slope of 2.62 and an outer slope of 3.8 at the break radius of 28 kpc, or an equally good Einasto profile (Einasto 1965) with a concentration index of 2.2 and effective radius of 22.2 kpc. Deason et al. (2011) analyzed ∼20,000 A-type photometric stars selected from the SDSS data release 8 (Ahn et al. 2012) and obtained a best-fitting BPL density with an inner slope of 2.3 and an outer slope of 4.6, with a break radius at 27 kpc and a constant flattening of 0.6. Subsequently, Deason et al. (2014) found a very steep outer halo profile with a power law of r−6 beyond 50 kpc, and yet steeper slopes of α = 6–10 at larger radii.
In addition, several pieces of work point to variations of the stellar halo flattening with radius. Preston et al. (1991) found that the density distribution of RR Lyrae follows a power law with α ∼ 3.2, together with a variable flattening changing linearly from 0.54 at center to 1 at 20 kpc. Subsequent work (Deason et al. 2011; Sesar et al. 2011) found evidence for flattening, but no evidence for a change with radius.
Spectroscopic maps of the stellar halo beyond ∼20 kpc in practice require luminous post-MS stars, as turnoff or other MS stars are too faint for current wide-field spectrographs. RR Lyrae and BHB stars have repeatedly been used as tracers to study the halo density profile, because they have precise distances and are bright enough to be observed at radii out to ∼100 kpc (Xue et al. 2008, 2011; Sesar et al. 2010; Deason et al. 2011, 2014). Yet, such stars, most prevalent in particularly old and metal-poor populations (Bell et al. 2010; Xue et al. 2011), are known to have a different structure and profile from the metal-poor red giant branch (RGB) stars (K giants). Moreover, Preston et al. (1991) and Morrison et al. (2009) claimed that BHB stars in the inner halo had a different radial profile from RR Lyraes. To this end, it is crucial to construct the halo shape and radial profile of the stellar halo in K giants. Xue et al. (2014) presented a catalog of K giants with unbiased distance estimates good to 16%, metallicities, velocities, and photometric information, drawn from the Sloan Extension for Galactic Understanding and Exploration (SEGUE; Yanny et al. 2009), which contains ∼150 stars beyond 60 kpc. SEGUE was focused on the halo in its targeting and has two phases, SEGUE-1 and SEGUE-2. The target selection of K giants changed a lot during SEGUE-1, but SEGUE-2 adopted a unified K-giant selection, so one can understand and model their selection function well for K giants observed in SEGUE-2. Owing to the well-understood selection function, it is possible to determine the halo profile and shape of these tracers, which is the main goal of the present paper. Specifically, we set out to describe the stellar halo distribution, presuming that the density is stratified on (oblate) spheroids, with a radial profile from 10 to 80 kpc that can be characterized by simple functional forms (either an Einasto profile or BPL). We also explore the metallicity dependence of the shape and radial profile of the stellar halo.
In the next section, we lay out the properties and the selection function of the SEGUE K giants. In Section 3, we present the method of fitting a series of parameterized models to SEGUE-2 K giants, explicitly and rigorously considering the selection function. This step is key in obtaining accurate radial profiles and metallicity gradient model in the meantime. The results for the stellar halo's radial profile flattening are presented in Section 4, along with the metallicity gradient in the stellar halo. Finally, Section 5 discusses the comparison between our results and previous work, as well as implications for dynamical models.
2. SEGUE-2 K GIANTS AND THEIR SELECTION FUNCTION
The SDSS (York et al. 2000) is an imaging and spectroscopic survey covering roughly a quarter of the sky, which has both ugriz imaging (Fukugita et al. 1996; Gunn et al. 1998; Stoughton et al. 2002; Pier et al. 2003; Eisenstein et al. 2011) and low-resolution spectra (λ/Δλ ∼ 2000). SEGUE is one of the key projects and has two phases, SEGUE-1 and SEGUE-2, which aim to explore the nature of stellar populations from 0.5 to 100 kpc (Yanny et al. 2009; C. M. Rockosi et al. 2015, in preparation). SEGUE-2 spectroscopically observed approximately 120,000 stars, focusing on the stellar halo of the Galaxy.
To understand the underlying spatial distribution of the K giants on the basis of this spatially incomplete sample, we need to understand (and account for) the probability that a star of a given luminosity, color, and metallicity ends up in the sample, given its direction and distance. Spectroscopic surveys of the Milky Way are inevitably affected by such selection effects (see Rix & Bovy 2013), often referred to as "selection biases." They arise from a set of objective and repeatable decisions of what to observe, necessitated by the survey design. In particular, only a small fraction of the sky was covered by SEGUE plates, and for most plates only a fraction of stars that satisfy the photometric selection criteria could be targeted with fibers. Finally, not all targeted stars yield spectra good enough to result in a catalog entry, i.e., had signal-to-noise ratio high enough to verify that they are giants and yield a metallicity measurement. Bovy et al. (2012) and Rix & Bovy (2013) spelled out how to incorporate this selection function in fitting a parameterized model for the stellar density, and we follow their approach in this and the next section.
Both the SEGUE-1 and SEGUE-2 surveys targeted halo giant star candidates, using a variety of photometric and proper-motion cuts. About 90% of the final K-giant sample came from objects observed as l-color K-giant targets. The l color (Lenz et al. 1998) is a photometric metallicity indicator for stars in the color range , designed to select metal-poor K giants.5 As mentioned in Section 1, only the selection function for K giants observed in SEGUE-2 can be understood well because SEGUE-2 adopted a consistent color–magnitude cut to select K giants throughout its entire survey: 15.5 < g0 < 18.5, r0 > 15, , , , and . Therefore, we restrict our analysis to this category. We also require that l-color K-giant candidates have good proper-motion measurements ≤11 mas yr−1. Since broadband photometry is a poor MS versus giant discriminator, not all stars targeted under the above criteria will be giants. The subsequent identification of K giants is based solely on their spectroscopic properties. As described in Xue et al. (2014), this requires spectra that have good Mg index and stellar atmospheric parameters determined by the SEGUE Stellar Parameter Pipeline (Lee et al. 2008a, 2008b, 2011) but have no strong G band.
To understand the selection function, we must compare the color–magnitude distribution of these spectroscopically confirmed K giants with the analogous distribution of all possible photometric l-color K giant candidates. Figure 1 shows these two distributions (as contours and grayscale, respectively), summed over all SEGUE-2 plates. As the marginalized histograms at the sides of the panels show, these two distributions are nearly indistinguishable: the chance of a photometric candidate being confirmed as a K giant is independent of color and magnitude. This simplifies the subsequent analysis and is testament to SEGUE-2's consistency of target selection, targeting, and spectral analysis.
However, while the selection function is constant with apparent magnitudes and dereddened colors, it varies from plate to plate, in particular with the Galactic latitude of the plate. Given the pencil-beam nature of the SEGUE survey, it makes sense to specify the selection fraction plate by plate. For each plate we define the number of spectroscopically confirmed K giant as Nspec, and the number of l-color K-giant candidates in the plate (both those that were targeted and those that were not) as Nphot. Thus, the plate-dependent selection function (shown as Figure 2) is given by
Download figure:
Standard image High-resolution imageAs we want to analyze the spatial distribution, we also restrict our sample to SEGUE-2 l-color K giants that have reliable distance estimates from Xue et al. (2014). To eliminate the contamination from the disk component, we cull K giants with [Fe/H] > −1.2 and kpc, which leads to a final sample of 2413 l-color K giants. Figure 3 illustrates the basic sample properties: its sky coverage and spatial distribution (without accounting for the selection function), and the distribution of metallicity against distance from Xue et al. (2014). The stars' distribution reflects the pencil-beam pattern of the SEGUE survey; their galactocentric distances range from 7 to 85 kpc; their mean metallicity is −1.75 dex, with some being as metal-poor as −3.5. The reader should note that these metallicities are not typical of the halo as a whole because of our choice to exclude all stars with [Fe/H] > −1.2.
Download figure:
Standard image High-resolution image3. MODELING THE STELLAR DENSITY AND METALLICITY DISTRIBUTION IN THE HALO
In this section we lay out a forward-modeling approach to describe the spatial and metallicity distribution of the stellar halo with a set of flexible but ultimately smooth and symmetric functions. Both the modeling practicalities and the astrophysics require that we fit the spatial distribution and the abundance distribution simultaneously. We defer to Section 4 the question of how sensitively such a "smooth" description depends on the question of including or excluding stars that are presumably members of recognizable substructure (see Belokurov et al. 2006; Bell et al. 2008).
3.1. A Parameterized Model for the Observables
We presume that the stellar halo distribution can be sensibly approximated by a spheroidal distribution with a parameterized radial profile, allowing for radial variations in the metallicity distribution, . The combination of pH and pFe denotes the halo parameters, and is the basic galactocentric radial coordinate, given the flattening q(r).
In this section, we spell out a straightforward and rigorous approach to determine the posterior probability distributions for halo parameters in light of the above data, our knowledge of the SEGUE-2 selection function, and well-established astrophysical priors on the luminosity function of giant stars. The number of halo parameters depends on the complexity of the model stellar halo distribution. This approach essentially follows Bovy et al. (2012) and Rix & Bovy (2013).
Since we already have good estimates of the distance modulus and the r-band absolute magnitude Mr for all objects in the sample (Xue et al. 2014), we treat as the observables defining , rather than , as it makes the fitting formalism more intuitive. However, the r-band apparent magnitude m and (dereddened) color c (referring in particular to ) will appear explicitly in the observational selection function. We denote the angular selection function as (see Equation (1)), the magnitude–color selection function by , the metallicity selection function by S([Fe/H]), additional spatial cuts by , and the luminosity selection function as , expressed in terms of the "observables" above. We denote the prior external information on the distribution of absolute magnitude on the giant branch as p(Mr); for old and metal-poor populations this is well established through cluster luminosity functions (see, e.g., Xue et al. 2014). The luminosity selection function will appear explicitly, because stars low on the giant branch were removed (depending on [Fe/H] ) in Xue et al. (2014) to avoid confusion between RGB and red clump (RC) stars. Figure 4 shows the limits of Mr for a given [Fe/H]. In short, we need to spell out and then quantify all terms that matter for predicting the rate function, i.e., the expected frequency or probability of finding a sample member with a given if a halo with pH and pFe were true.
Download figure:
Standard image High-resolution imageGiven pH and pFe, the expected rate function for finding a star with is then
This rate function is simply a concise way of specifying what set of observations we expect, given an assumed halo model, including all observational selection effects and pertinent prior information. The price for putting this in one place is that the rate function has quite a number of distinct terms. In the next two subsections, we spell out and discuss these terms. To make the rate function a probability, it must be normalized for every new set of trial model parameters (pH, pFe); this is the only time-consuming step in such modeling.
3.2. Models for the Spatial and [Fe/H] Distribution of the Stellar Halo
Following a number of previous studies, we presume that the overall radial density profile of the halo can be described by an Einasto profile or by a (multiply) broken power law, with the density stratified on surfaces of constant rq in all cases; in addition, we devise a simple model for radial [Fe/H] variations.
3.2.1. Einasto Profile
The Einasto profile (Einasto 1965) is the 3D analog to the Sérsic profile (Sérsic 1963) for surface brightnesses and has been used to describe the halo density distribution (Merritt et al. 2006; Deason et al. 2011; Sesar et al. 2011):
where ν0 is the (here irrelevant) normalization, reff is the effective radius, n is the concentration index, and , for n ≥ 0.5. The free parameters of an Einasto profile with a constant flattening are .
3.2.2. Broken Power-law Profiles
Broken (or even multiply broken) power laws (BPLs) are another family of functional forms that has been used extensively to described the radial profile of the Galactic stellar halo (Saha 1985; Deason et al. 2011, 2014; Sesar et al. 2011). In most cases the change in the power-law index, , has been taken to be a step function. We adopt
In addition to the flattening and the (irrelevant) normalization, a (singly) broken power law has three parameters: αin, αout, and rbreak. Of course, this profile family encompasses an SPL, and it can be generalized to a multiply broken power law (twice-broken power law (TPL)) by introducing an additional pair of (see Deason et al. 2014).
3.2.3. Halo Flattening
While Preston et al. (1991) found evidence for a decrease of flattening with increasing radius, others did not (Deason et al. 2011; Sesar et al. 2011). As there is evidence that at least the innermost part of the halo is quite flattened (Carollo et al. 2010), we explore how our sample can inform us about radial variations in the flattening of the stellar halo beyond rGC = 10 kpc. To date no particular functional form to parameterize the possible variation of halo flattening has been established; therefore, we consider the functional form for q(r) as
where q0 is the flattening at center, changing to at large radii, with r0 the (here exponential) scale radius, over which the change of flattening occurs.
3.2.4. Radial Variations in the Metallicity Distribution
The existence of a halo metallicity gradient is currently controversial (see, e.g., Carollo et al. 2007; Schönrich et al. 2011; Fernández-Alvar et al. 2015). This is partially due to the lack of in situ tracers for which accurate abundance measures are possible with low-resolution spectroscopy, and partially due to the lack of attention to the effects of the complex SEGUE selection function. This is the first paper that has applied a forward-modeling approach to the K giants in the halo, addressing both of these concerns. We choose a mixture model for the halo metallicity distribution as follows. The metal-rich component can be described by a Gaussian with a mean at −1. 4 dex and a sigma of 0.2 dex, and the metal-poor component follows a Gaussian with a mean at −2.1 dex and a sigma of 0.35 dex. Therefore, we suppose the metallicity distribution, , as a radially varying combination of these two metallicity components. The metallicity distribution model is also expressed as a function of rq. We choose
where rq is the same as defined in Equation (3); (mean, dispersion) are Gaussian functions. Besides the common parameter q, the metallicity gradient has two other parameters, f0 and γ (we mark them as pFe).
3.3. Selection Effects and Prior Information
The remaining terms in Equation (2) incorporate the various selection effects and pieces of prior information into the prediction of the rate function and are given in the following.
The Jacobian term is given by
where Ωplate = 7 deg2 and we presume the Jacobian to be constant within a plate. The prior on the absolute magnitude distribution along the giant branch (Xue et al. 2014) is given by
The sample's metallicity cuts, aimed at eliminating bulge and thick-disk contributions, as well as any spurious metallicity determinations, −3.5 < [Fe/H] < −1.2, are reflected in
The spatial cuts to geometrically excise any bulge and thick-disk stars are
The reliable distance determination requires minimizing the contamination of RGB stars by RC stars, which translates into [Fe/H]-dependent absolute magnitude cuts (Figure 4):
The set of observed SEGUE plates leads to a rather sparse angular selection function in (l, b), illustrated in Figures 2 and 3:
Finally, the apparent magnitude and color cuts (see Figure 1) for the SEGUE targeting are reflected in
Taking these terms together and following Bovy et al. (2012), we can now directly calculate the likelihood of the data given (pH,pFe) and the rate function:
where i is the index of each K giant and we use that the data points are identically and independently drawn. The normalization cλ is the integral over the volume in the space,
Here we have already performed the dldb integral for each plate. This normalization integral is the computationally most expensive part of the parameter estimates. But it can be computed efficiently using Gaussian quadratures, where we adopt 48 × 20 × 20 transformation points in , and [Fe/H] space, and where the parameter-independent parts such as the Jacobian term, the luminosity prior, and selection functions can be pre-computed on a dense grid. We assume the priors on all parameters (pH, pFe) to be flat over the pertinent range; therefore, the posterior distribution function (pdf) of the parameters, , is proportionate to the likelihood (Equation (14)), differing only in its units ("1/parameters" versus "1/data"). We then vary the (pH, pFe) to sample the parameter pdf using emcee, which implements an efficient Markov chain Monte Carlo technique6 (Foreman-Mackey et al. 2013).
4. RESULTS
We now present the results of applying the modeling from Section 3 to the sample of Section 2. We illustrate these results in two ways: first, by showing the joint pdf's of the halo model parameters; second, we show what the radial and flattening profiles actually look like, by drawing samples from these parameter pdf's. We start out with the simplest model (the one with the fewest parameters), an Einasto profile of constant flattening, to illustrate the results, present the basic gradient in the [Fe/H] distribution, and explore the impact of excising stars that are manifestly in substructures. We then proceed to include variable flattening and the BPLs as radial profiles. The results are summarized in Table 1.
Table 1. A Summary of Our Best-fitting Models
Model | Parameters | Np | a | ΔBICb |
---|---|---|---|---|
Excl. all substructures | ||||
[Fe/H]-model | f0 = 0.6 ± 0.01, γ = −0.2 ± 0.05 | 2 | ⋯ | ⋯ |
+ | ||||
SPL | α = 3.6 ± 0.1, q = 0.68 ± 0.02 | 2 | −35 | 54 |
Einasto | n = 3.1 ± 0.5, reff = 15 ± 2 kpc | 3 | −16 | 25 |
q = 0.7 ± 0.02 | ||||
BPL | αin = 2.1 ± 0.3, αout = 3.8 ± 0.1 | 4 | −14 | 27 |
rbreak = 18 ± 1 kpc, q = 0.7 ± 0.02 | ||||
TPL | α1 = 2.1 ± 0.2, α2 = 3.8 ± 0.1 | 4 | −13 | 26 |
α3 = 4.8 ± 0.8, q = 0.7 ± 0.02 | ||||
rbreak1 = 18 kpc, rbreak2 = 65 kpc | ||||
[Fe/H]-model | f0 = 0.6 ± 0.01, γ = −0.26 ± 0.06 | 2 | ⋯ | ⋯ |
+ | ||||
Einasto-q(r) | n = 7.7 ± 1.5, reff = 2 ± 1 kpc | 5 | −1 | 9 |
q0 = 0.2 ± 0.1, qinf = 0.78 ± 0.05 | ||||
r0 = 6 ± 1 kpc | ||||
BPL-q(r) | αin = 4.2 ± 0.4, αout = 4.3 ± 0.3 | 6 | 0 | 14 |
rbreak = 22 ± 23 kpc, q0 = 0.2 ± 0.1 | ||||
qinf = 0.8 ± 0.03, r0 = 6 ± 1 kpc | ||||
SPL-q(r) | α = 4.2 ± 0.1, q0 = 0.2 ± 0.1 | 4 | 0 (−12758) | 0(25546) |
qinf = 0.8 ± 0.03, r0 = 6 ± 1 kpc | ||||
Incl. all substructures | ||||
[Fe/H]-model | f0 = 0.6 ± 0.01, γ = −0.15 ± 0.04 | 2 | ⋯ | ⋯ |
+ | ||||
SPL | α = 3.4 ± 0.1, q = 0.74 ± 0.02 | 2 | −62 | 110 |
Einasto | n = 2.3 ± 0.2, reff = 18 ± 1 kpc | 3 | −24 | 41 |
q = 0.77 ± 0.02 | ||||
BPL | αin = 2.8 ± 0.1, αout = 4.3 ± 0.1 | 4 | −25 | 52 |
rbreak = 29 ± 2 kpc, q = 0.77 ± 0.02 | ||||
TPL | α1 = 2.8 ± 0.1, α2 = 4.3 ± 0.2 | 4 | −26 | 53 |
α3 = 4.2 ± 0.5, q = 0.77 ± 0.02 | ||||
rbreak1 = 30 kpc, rbreak2 = 55 kpc | ||||
[Fe/H]-model | f0 = 0.62 ± 0.01, γ = −0.22 ± 0.06 | 2 | ⋯ | ⋯ |
+ | ||||
Einasto-q(r) | n = 6.1 ± 1.7, reff = 3 ± 2 kpc | 5 | 0 (−17675) | 6 |
q0 = 0.3 ± 0.05, qinf = 0.9 ± 0.04 | ||||
r0 = 9 ± 2 kpc | ||||
BPL-q(r) | αin = 4.2 ± 0.3, αout = 4.5 ± 0.4 | 6 | 0 | 12 |
rbreak = 32 ± 18 kpc, q0 = 0.3 ± 0.1 | ||||
qinf = 0.9 ± 0.04, r0 = 9 ± 1 kpc | ||||
SPL-q(r) | α = 4.4 ± 0.1,q0 = 0.3 ± 0.1 | 4 | −2 | 0(35384) |
qinf = 0.9 ± 0.04, r0 = 9 ± 1 kpc |
Notes. We give the type of model, the best-fitting parameters of the model, the number of free parameters, difference in log-likelihood from maximum likelihood value, and difference in BIC from minimum BIC value. Parameters that are kept fixed are highlighted in bold.
a = -, so means that the model has maximum likelihood. The value in parentheses is . b , and ΔBIC = , so the best model has ΔBIC = 0. The value in parentheses is .Download table as: ASCIITypeset image
Figure 5 illustrates the result of the simplest model we fit to the entire data set, by sampling the parameter pdf using Equation (14). This model has three parameters for , , and two parameters for the metallicity distribution , namely . The projections in Figure 5 of the pdf involving only pH are shown in gray, while those also involving pFe are shown in blue. This figure illustrated that the sample size and data quality of the sample are sufficient to provide formally very well constrained parameters.
Download figure:
Standard image High-resolution image4.1. Radial Gradients in the Halo's Metallicity Distribution
Figure 5 shows that the mix of the intermediate () and metal-poor () [Fe/H] components changes with radius (i.e., ): the relative importance of the metal-poor component increases toward large radii, by a factor of 1.4 ± 0.1 from 10 to 65 kpc. Our best-fit metallicity distribution predicts a mean metallicity of −1.8 in this radial range. The figure also shows that the covariances between the fits to the stellar density and to the [Fe/H] distribution are weak. The resulting [Fe/H] distribution, projected into the space of the rate function, is shown in the inset of Figure 5; this inset illustrates that the distribution of the actual sample members in the plane is plausible. Note that nearly all sample stars within rGC = 10 kpc, much of the "inner halo" of Carollo et al. (2007), have been eliminated from the fit (see Figure 3); therefore, our result is the detection of an outward metallicity gradient within the "outer" halo.
4.2. The Impact of Halo Substructures on Fitting Smooth Models
As discussed in Section 1, both cosmological models and observations imply that a good portion of halo stars, at least beyond 20 kpc, are in substructures. Especially the prominent ones, such as the Sagittarius stream and the Virgo overdensity, can and will affect the fits of smooth models, as pointed out by Deason et al. (2011). Recently, Janesh et al. (2015) used a position–velocity clustering estimator (the 4-distance) in combination with a friends-of-friends algorithm to identify the stars in the Xue et al. (2014) K-giant sample belonging to substructures; they found that 27% of the K giants are clearly associated with the Sagittarius streams, the Orphan streams, the Cetus Polar stream, and other unknown substructures. The results of Janesh et al. (2015) provide a straightforward algorithmic way of excising much of the manifest substructure from the sample. For the Sagittarius stream (Belokurov et al. 2014) we know the expected position–velocity distribution quite well, which suggests that six of the most distant sample members are likely members of Sagittarius's trailing arm. These lie at rGC > 60 kpc, , and cluster tightly around Vgsr = 100 kms−1, which also eliminates them from the original sample of 2413 l-color K giants, leaving 1757 l-color K giants.
We then repeat the fitting of the same model as illustrated in Figure 5 and find—unsurprisingly—significant differences (Figure 6): an Einasto profile that is more concentrated, has a smaller effective radius, and is more flattened, n = 3.1 ± 0.5, reff = 15 ± 2 kpc, and q = 0.7 ± 0.02. The metallicity gradient, however, remains very similar with f0 = 0.6 ± 0.01 and γ = −0.2 ± 0.04. Such a model already provides quite a good fit to the data, as Figure 7 shows: the model prediction for the distribution of the sample members, averaged over all directions and metallicities and accounting for all selection effects, matches the actual distribution quite well. This figure corresponds to the top right inset in Figure 6, just marginalized over metallicities.
Download figure:
Standard image High-resolution imageDownload figure:
Standard image High-resolution image4.3. Different Functional Forms for the Radial Profile
Following previous approaches, we next explore different functional forms for the radial profile, assuming radially constant flattening.
Figure 8 shows the results analogous to Figure 6, but for a BPL: at a best-fit, constant flattening of q = 0.7 ± 0.02 we find an inner slope of αin = 2.1 ± 0.2, distinctly different from the outer slope of αout = 3.8 ± 0.1, with a break radius of rbreak = 18 ± 1 kpc. The flattening q is consistent with that of the best-fit Einasto profile. Again, also with a BPL profile remains similar and is consistent with the observations. The Bayesian information criterion (BIC; Schwarz 1978; also called the Schwarz criterion) is a criterion for model selection. It is defined as , where is the maximized likelihood of the model, Np is the number of parameters in the model, and Ndata is the sample size. BIC takes into account both the statistical goodness of fit and the number of parameters that have to be estimated to achieve this particular degree of fit, by imposing a penalty for increasing the number of parameters. The model with the lowest BIC is preferred. As the values of ΔBIC in Table 1 indicate, the Einasto and BPL radial profiles fit the data comparably well. A TPL, with the break radii held fixed at the radii suggested by Deason et al. (2014), leads to a comparably good fit (see Figure 9 and Table 1). Compared to Deason et al. (2014), the fit is consistent within 65 kpc; it shows no evidence for a steep drop beyond, but our present sample is vary sparse at these distances. For reference, we also fit the data with an SPL of constant flattening, but on the basis of its ΔBIC it can be ruled out compared to the BPL and TPL models. Our results confirm that if constant flattening is assumed, the halo profile must be described by a radial profile break at ∼20 kpc.
Download figure:
Standard image High-resolution imageDownload figure:
Standard image High-resolution image4.4. Radial Variation of Halo Flattening
We now proceed to allow more model flexibility by allowing the flattening q(r) to vary with radius, according to Equation (5); we fit both the Einasto profile and power laws (BPL and SPL) to the data. Allowing for flattening variations makes for distinctly better fits to the data, as quantified by ΔBIC (Table 1). These fits imply a strong variation of the halo flattening with radius, as the pdf's for an Einasto profile in Figure 10 illustrate: while at large radii , the flattening at asymptotically small radii becomes formally very strong, q0 ≈ 0.2. Given that the characteristic radius at which this flattening change occurs is kpc, and hence smaller than the minimal radius of all data points at 10 kpc, the actual flattening changed across the radial range constrained by the data is more modest, as Figure 12 illustrates: q changes from 0.55 ± 0.02 at 10 kpc to 0.8 ± 0.03 at large radii. Allowing for variable flattening also changes the Einasto parameters quite drastically: the formal effective radius becomes very small, and there is a strong covariance between and concentration n. However, the actual predictions for the slope of the radial profile, , in the radial regime constrained by the data are quite similar between models with constant and variable flattening, as Figures 13 and 14 show. This just illustrates to compare model fits drawn from different functional form families in their projection in the space of observables, not just in the space of the parameters. The analogous pdf(s) of the analogous fit for a BPL look at first sight bewildering. However, it reveals that the implied flattening profile is the same as for an Einasto fit. The unconstrained pdf's on the break radius simply reflect the fact that the inner and outer power-law slopes become indistinguishable, which is equivalent to an SPL. Therefore, we turn to fit the data to an SPL with a varying flattening (see Figure 11). Remarkably, we find a very good fit with an SPL in rq when allowing for variable flattening, because this model has the minimum BIC, as illustrated by ΔBIC = 0 in Table 1. The flattening profile is the same as for an Einasto fit, as also shown in Figure 12: there is effectively a break in the flattening profile at 15–20 kpc. Figure 13 shows that the actual local density slope along an intermediate axis for this SPL is similar to the other fits. Taken together, this implies in the context of these density models that the stellar halo density is best and simplest described (10–80 kpc) by an SPL in the variable rq, with the break in the radial profile at 20 kpc replaced by a break in the flattening at a comparable radius (see Figure 12).
Download figure:
Standard image High-resolution imageDownload figure:
Standard image High-resolution imageDownload figure:
Standard image High-resolution imageDownload figure:
Standard image High-resolution image5. DISCUSSION AND SUMMARY
Using a large spectroscopically confirmed sample of K giants, well-understood population tracers of the stellar halo, we have attempted to characterize the radial profile, the flattening profile, and the metallicity profile of the "outer halo" from 10 to 80 kpc; we have done so after excising algorithmically identified members of distinct halo substructures from the sample. On the one hand, this analysis has re-emphasized that the choice of the functional forms for the density models matters in casting the results; on the other hand, we could show that three robust aspects emerge: the profile of the local density slope, , when taken along an axis intermediate between major and minor axis, is consistent among all functional profile forms we explored; the stellar halo distribution appears distinctly flatter within 20 kpc; and there is a slight, but significant, outward metallicity decrease within the "outer halo." We also find a break in the properties of the halo density distribution at ∼20 kpc. This can be attributed to a break in the radial power law, as in a number of previous analyses, but we find that it foremost reflects a break in the halo flattening: the data are well fit by an SPL in the spheroidal coordinate rq, with a distinct change in q(r) at a radius of ∼20 kpc.
We illustrate the first point about the local slope in Figure 13. As mentioned in the Introduction, this slope plays a particular role in dynamical modeling, perhaps most easily seen in the case of axisymmetric Jeans equation modeling. Figure 13 compares the radial profile functions between the models with constant or variable flattening: all models have consistent radial profile slopes, , within 65 kpc. There is remarkable qualitative agreement, in the sense that these slopes vary consistently with radius, within the constraints of the functional form imposed. In all cases does the halo profile steepen between 10 and 65 kpc. Clearly, the fits with fixed q differ more from one another and from the q(r) fits, as they must through their more stringent functional form constraints. Allowing for a yet steeper drop beyond 65 kpc with a TPL profile leads to indistinguishable results. The largest discrepancies are between power laws and Einasto profiles beyond 65 kpc; we attribute that to the Einasto profile being constrained by the abundance of data at rgc ≤ 65 kpc (see bottom panel of Figure 13), which inevitably leads to at large radii. A look at Table 1 shows, however, that all of these models make the data appear comparably likely; the likelihood differences are significant, but not drastic (see also Figure 7). In particular, the Einasto and BPL profiles lead to near-identical data likelihoods; the model fits with the flattening forced to be constant make the data significantly less likely. Remarkably, the SPL with variable flattening has nearly the same along the intermediate axis, where the seeming radial change in power-law index is simply a reflection of the change in the radial coordinate rq, attributable to q(r).
We also compare these findings to other work, in particular that of Deason et al. (2011, 2014). Figure 14 shows that the halo profile, as traced in this analysis via K giants, is consistent with the findings of Deason et al. (2011) within 65 kpc. To follow up the findings of Deason et al. (2014), a steep drop in the density of BHB stars beyond 65 kpc, we specifically fit a TPL with parameters (α1, α2, α3, q). For comparison, we held the break radii fixed at rbreak1 = 18 kpc and rbreak2 = 65 kpc. Yet, as Figure 14 reveals, the halo slopes in BHB stars and found here in K giants are formally inconsistent beyond 65 kpc. Yet, there are only seven stars beyond 65 kpc, so the paucity of distant K giants precludes a more stringent comparison. Besides, another two avenues are conceivable: first, we know that the substructure differs between BHB stars and giants and MS turnoff stars (Bell et al. 2008; Xue et al. 2011; Janesh et al. 2015). Second and related, it is conceivable than many of the stars at large radii in the present K-giant sample could be Sagittarius stream members, even though we made an attempt to remove such stars from the sample.
Download figure:
Standard image High-resolution imageIn summary, we believe that the present study has brought forward a number of new aspects: First, we have worked out a forward modeling of the spectroscopic data that has not been applied previously in this context. We believe that, in particular for tracers that are not standard candles like BHB stars, such an approach is warranted and powerful. Second, we were able to show on this basis (1) that there is an outward metallicity gradient in the halo beyond 10 kpc; (2) that the halo is distinctly flatter between 10 and 20 kpc, compared to larger radii; this distinct change in flattening suggests that it is more appropriate to think of the break in halo profile at 20 kpc as a break in flattening, rather than as a break in the radial profile at forced constant flattening; and (3) that there is overall consistency with previous analyses when it comes to the dynamically important quantity in the radial range 10–65 kpc.
The research has received funding from the European Research Council under the European Union's Seventh Framework Programme (FP 7) ERC grant agreement No. [321035]. X.-X.X. acknowledges the Alexander von Humboldt Foundation for a fellowship that enabled this work and the support from NSFC grants 11103031, 11233004, 11390371, and 11003017. H.W.R. acknowledges funding by the Sonderforschungsbereich SFB 881 The Milky Way System (subproject A3) of the German Research Foundation (DFG). J.B. acknowledges support from a John N. Bahcall Fellowship and the W. M. Keck Foundation. H.L.M. acknowledges support from NSF grant AST-1009886. We thank Glenn van de Ven and Ling Zhu for helpful discussions, and Vasily Belokurov and Alis Deason for valuable input on an earlier version of the manuscript. Funding for SDSS-III has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, and the U.S. Department of Energy Office of Science. The SDSS-III Web site is http://www.sdss3.org/. SDSS-III is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS-III Collaboration, including the University of Arizona, the Brazilian Participation Group, Brookhaven National Laboratory, University of Cambridge, Carnegie Mellon University, University of Florida, the French Participation Group, the German Participation Group, Harvard University, the Instituto de Astrofisica de Canarias, the Michigan State/Notre Dame/JINA Participation Group, Johns Hopkins University, Lawrence Berkeley National Laboratory, Max Planck Institute for Astrophysics, Max Planck Institute for Extraterrestrial Physics, New Mexico State University, New York University, Ohio State University, Pennsylvania State University, University of Portsmouth, Princeton University, the Spanish Participation Group, University of Tokyo, University of Utah, Vanderbilt University, University of Virginia, University of Washington, and Yale University.
Footnotes
- 5
- 6
The user guide for emcee can be found at http://dan.iel.fm/emcee/current/.