Chemically Dissected Rotation Curves of the Galactic Bulge from Main-sequence Proper Motions^*

In addition to velocity dispersion trends, the trend in bulge mean radial velocity (against Galactic longitude or galactocentric radius) might also be expected to vary with metallicity, but here the magnitude (or even existence) of such a dependence is less clear. Earlier spectroscopic surveys suggest a clear difference between metal-poor and metal-rich samples. For example, Harding & Morrison (1993) and Minniti (1996) demonstrated that "metal-rich" stars show a gradient in circular speed with galactocentric radius, consistent with the "solid-body"-type rotation traced by planetary nebulae (Kinman et al. 1988), Mira variables (Menzies 1990), and SiO masers (Nakada et al. 1993). In contrast, metal-poor objects (using [Fe/H] ≲ −1.0, and thus likely including a large contribution from the inner halo) showed no strong evidence for a rotational trend. More recently, Kunder et al. (2016) found that their metal-poor RR Lyrae sample with mostly subsolar metallicities (−2.4 < [Fe/H] ≲ +0.3, peaking at [Fe/H] ∼ −1.0) shows no strong signature of rotation from radial velocities in any Galactic latitude range. This is in contrast to the majority-bulge population, which shows bulk rotation with amplitude v_GC ± ≈80 km s⁻¹ progressing from the first to fourth Galactic quadrant (e.g., Howard et al. 2009; Kunder et al. 2016). This rotation-free component is estimated to be a rather small part of the overall bulge stellar population (Kunder et al. 2016).

Restricting attention to [Fe/H] ≳ −1.0 (to sample mainly bulge and disk stars), the body of more recent spectroscopic studies does not show strong evidence for metallicity dependence of radial velocity rotation curves (usually plotted against Galactic longitude). For example, the ARGOS survey (Ness et al. 2013b) and the Gaia-ESO survey (Williams et al. 2016) each show no strong difference between metal-rich and metal-poor bulge objects (the studies use slightly different cuts for metal-rich and metal-poor objects). However, the Giraffe Inner-Bulge Survey (GIBS, which is unusual among the spectroscopic studies in reaching as close as b = −2° to the Galactic midplane) shows a possible difference in rotation curve slope between objects at $[\mathrm{Fe}/{\rm{H}}]\lt -0.3$ and [Fe/H] > +0.2; however, at about 1.5σ significance, the difference is not yet compelling (Zoccali et al. 2017). Thus, the radial velocity surveys focusing on the majority-bulge population (with [Fe/H] ≳ −1.0) show no strong metallicity dependence in the trends of mean radial velocity against Galactic longitude.

Proper motions offer an independent method to kinematically chart the bulge rotation curves and, if information on chemical composition is available, explore whether multiple abundance samples really do show distinct mean motions, as well as the well-established velocity dispersion differences.

To date, proper-motion investigations in the context of multiple populations (or a continuum) have mostly been performed using bright giants. For example, in addition to the vertex deviation investigations reported in the previous section, proper motions of bulge giants using OGLE (Poleski et al. 2013) and with the Wide Field Imager on the La Silla 2.2 m telescope (Vásquez et al. 2013) have been used to uncover azimuthal streaming in the bulge X-shaped structure. However, Qin et al. (2015) caution via N-body models that the underlying bar pattern speed cannot directly be constrained just from the near-side/far-side longitudinal proper-motion difference.

The above radial velocity and proper-motion studies all use bright giants as tracers, often red clump giants (RCGs), which are much less spatially crowded from the ground than are main-sequence (MS) objects. This causes them to be limited by the small intrinsic population size per field of view. For example, ARGOS typically observed about 600 stars at [Fe/H] > −1.0 per 2^◦-diameter field of view (Ness et al. 2013b). Thus, mean velocities interpreted for rotation trends represent averages both over quite large angular regions on the sky and, more importantly, over the entire distance range along the line of sight.

To make further progress, an independent measure of bulge rotation is needed, using a tracer sample sufficiently populous that the sample can be dissected by line-of-sight distance to mitigate the statistical limitations of giant branch tracers. MS tracers are orders of magnitude more common on the sky, affording the opportunity to dissect a single sight line along the line of sight, thus offering a "pencil-beam" complement to the wide-field surveys that use the bright end of the CMD.⁹

It is the charting of the chemically dissected bulge rotation curve from MS proper motions that we report here. Because this is a relatively new technique, we briefly review the short literature in MS proper-motion bulge rotation curve determination before proceeding further.

Proper-motion-based rotation curves¹⁰ from MS bulge stars are relatively rare in the literature.¹¹ Kuijken & Rich (2002) were the first to demonstrate the approach for MS populations, for both the Baade and Sagittarius windows, presenting the HST/WFPC2-derived rotation and dispersion curves against photometric parallax (with photometric parallax determined as a linear combination of color and magnitude in order to remove the color–magnitude slope of the MS tracer population of interest). This demonstrated a clear sense of rotation, with the near side of the bulge showing positive mean longitudinal proper motion relative to the far side (a determination made before the much brighter RCGs were used to show bulge rotation from proper motions; Sumi et al. 2004). The proper-motion dispersion showed a slight increase in the most populous middle bins of photometric parallax (most strongly pronounced in the latitudinal proper-motion dispersion σ_b) for their Sagittarius window field. Kuijken (2004) presented an extension of this work to multiple fields across the bulge, including the use of three minor-axis fields to estimate the vertical gravitational acceleration along the Galactic minor axis.

Kozłowski et al. (2006) were able to demonstrate similar behavior to the Kuijken & Rich (2002) rotation curves in their analysis of proper motions in Baade's window. This was the only field for which a sufficiently large sample of sufficiently precisely measured MS stars could be measured from their large 35-field study (which used WFPC2 for early-epoch observations and ACS/HRC for late-epoch observations). While their dispersion curve is consistent with a flat distribution, the rotation trend in Galactic longitude was clearly observed. Kozłowski et al. (2006) may also have been the first to detect the weak trend in latitudinal proper motion μ_b due to solar reflex motion (see Vieira et al. 2007, for discussion of this effect, including its detection using sets of ground-based observations of bulge giants over a 21 yr time baseline). In any case, Kozłowski et al. (2006) were the first to detect the proper-motion correlation C_l,b at statistical significance from any population (using the RCGs that formed their main target population), using it to constrain the tilt angle of the bulge velocity ellipsoid. As they point out, detection of C_l,b (or equivalently the orientation angle ϕ_lb of the proper-motion ellipsoid) allows constraints to be placed on the orbit families for bulge populations, although the conversion from observation to physical constraint is not simple (e.g., Zhao et al. 1994; Häfner et al. 2000; Rattenbury et al. 2007).

Clarkson et al. (2008, hereafter Cl08) extended the rotation curve approach, using a much deeper data set with ACS/WFC toward the Sagittarius window, estimating photometric parallax directly with reference to a fiducial isochrone describing the average population in the CMD. Consistent with Kuijken & Rich (2002) and Kozłowski et al. (2006), this showed a clear sense of rotation in Galactic longitude, a clear detection of the latitudinal proper-motion trend from near side to far side, and a pronounced peak in the velocity dispersion of both coordinates (σ_l and σ_b) coincident with the most densely populated section of the photometric distance range of the sample. Cl08 converted proper motions to velocities, charting the run of the mean velocity (i.e., the rotation curves), the semiminor and semimajor axis lengths (i.e., the velocity dispersions), and the variation of the orientation ϕ_lb of the projected velocty ellipse with line-of-sight distance, and verified through simulation and comparison with the behavior of RCGs that indeed distance effects are observable in MS photometric parallax (though unlike RCG tracers, unresolved binaries blur somewhat the inferred distances for a given MS population).

More recently, in a careful study of three off-axis bulge fields using WFPC2 for early-epoch observations and ACS/WFC for the late epoch, Soto et al. (2014) were able to extract the rotation curve (and associated proper-motion dispersion curves) for a field farther from the midplane, at (l, b) = (+358, −717).¹² Soto et al. (2014) also computed the run of velocity ellipse orientation ϕ_lb with photometric distance, finding trends consistent with Cl08. The kinematics of MS objects at some distance from the plane were thus established to be broadly similar to those at the more central Baade and Sagittarius window fields.

Ground-based surveys are now starting to measure proper motions for MS bulge objects. For example, proper motions from the VVV survey have already been used to draw proper-motion rotation curves for both giant branch and upper MS populations (although the upper MS population shows much higher proper-motion scatter and substantially different selection effects compared to the giants; Smith et al. 2018).

The lack of metallicity information for MS populations has limited both the measurement accuracy and scientific applicability of MS proper-motion rotation curves. The [Fe/H] spread for bulge populations contributes a scatter of up to ∼1 mag on the MS (e.g., Haywood et al. 2016), competing with the photometric parallax signal due to the intrinsic distance distribution along the line of sight. While comparison with the behavior of RCGs suggests that indeed the rotation curve can be recovered, a lack of [Fe/H] information for the MS tracers contributes to substantial mixing in photometric parallax that can dilute the signature of underlying rotation (Clarkson et al. 2008). Conversely, charting bulge proper-motion rotation curves from samples partitioned by relative metallicity allows an independent probe of the chemical and dynamical correlations resulting from the complex formation and evolutionary processes at work in the inner Milky Way.

Until recently, no observational data set existed that would allow the proper-motion-based rotation curves to be charted for multiple spatially overlapping MS metallicity samples in the bulge, as the relevant tracer samples (a few magnitudes beneath the MS turnoff, and well clear of the subgiant and giant branches in the CMD) are far too faint and spatially crowded for objects to be chemically distinguished using current spectroscopic technology.

The situation changed with the WFC3 Bulge Treasury Survey (hereafter BTS; Brown et al. 2009), which used three-filter flux ratios to construct a "temperature" index $[t]$ (a function of F555W, F110W, and F160W magnitudes, similar to V, J, H) and a "metallicity" index $[m]$ (using F390W, F555W, and F814W magnitudes, similar to Washington-C, V, I), with scale factors chosen so that $[t]$ and $[m]$ are relatively insensitive to reddening. This allows stars to be chemically tagged in a relative sense by their location in $[m]$, $[t]$ space, down to much fainter limits and in regions of higher spatial density than currently allowed by spectroscopy. Brown et al. (2010) showed that indeed the wide bulge metallicity range can be traced photometrically by this method, setting $[t]$ and $[m]$ indices for tens of thousands of MS objects in each of the four observed bulge fields. Inverting the photometric indices then produced relative [Fe/H] distributions broadly similar to the spectroscopic indications from much brighter objects (e.g., Hill et al. 2011; Johnson et al. 2013). Computing these indices appropriately for objects near the bulge MS turnoff, Brown et al. (2010) found that the candidate exoplanet hosts of the SWEEPS field (Sahu et al. 2006) tend to pile up at the metal-rich end of the $[m]$ distribution as expected, suggesting that $[m]$ is indeed tracking metallicity. Exploitation of this unique data set to directly constrain the star formation history of the bulge is ongoing (see Gennaro et al. 2015, for an example of the techniques involved).

Here we combine the relative metallicity estimates from WFC3 BTS photometry with ultradeep proper motions using ACS/WFC, to construct the proper-motion-based rotation curves of candidate "metal-poor" and "metal-rich" MS samples, and examine whether and how the kinematics of the two samples differ from each other. Our work represents the first extension of chemodynamical studies of the bulge down to the MS.

This paper is organized as follows. The observational data sets are introduced in Section 2, with the techniques used to classify samples as "metal-poor" or "metal-rich" and to draw rotation curves described in Section 3. The rotation curves themselves are presented in Section 4. Section 5 discusses the implications of our results for the distribution of populations within both the bulge and proper-motion sample selection and discusses the impact of various systematic effects, with conclusions outlined in Section 6. Appendices A–I provide supporting information, including the full set of results in tabular form.

The SWEEPS data set used here consists of an extremely deep imaging campaign with a 9 yr time baseline using ACS/WFC in F606W and F814W (programs GO-9750, GO-12586, and GO-13057; PI K. C. Sahu). The observations and analysis techniques used to produce the proper motions and photometry used herein are described in some detail in previous papers (Sahu et al. 2006, hereafter Sa06; Cl08; Calamida et al. 2014, hereafter Ca14; Calamida et al. 2015, hereafter Ca15; Kains et al. 2017). Here we briefly describe the relevant characteristics for the present study.

Stellar positions in individual images were estimated using the distortion solution and effective point-spread function (PSF) methods developed by J. Anderson for HST and implemented for ACS/WFC in the img2xym.F routine (Anderson & King 2006) and associated utilities. This yields highly precise position measurements in a reference frame that is nearly free of distortion. With these techniques, per-measurement random uncertainties are as small as ${\epsilon }_{X,Y}\approx 0.002$ pixels per coordinate (e.g., Figure 3 of Cl08), and residual distortion is as low as ∼0.01 pixels (Anderson & King 2006; see also Appendix A). Detailed discussion of the methods can be found in Anderson & King (2006) and Anderson et al. (2008a, 2008b).

The 2011–2012–2013 epoch consists of 60 (61) images in F606W (F814W) taken with an approximately 2-week cadence, while the 2004 epoch consists of 254 (265) exposures in F606W(F814W) taken over a 1-week interval in 2004 (Sa06, all exposures in both programs being ≈5.5 minutes each, which well samples the bulge MS and minimizes downtime for buffer dumps).

Because the disk and bulge stars move relative to each other, the 2011–2012–2013 images were reduced separately from those in the 2004 epoch. Proper motions were derived from the best-fit positional differences between the 2004 and 2011–2012–2013 data sets; they thus represent two-epoch proper motions, but with positions in each individual epoch measured to very high accuracy. The positional differences (in ACS/WFC pixels) were rotated into a frame aligned with Galactic coordinates and converted from a displacement in pixels into rate of positional change in mas yr⁻¹ using the ACS/WFC plate scale (50 mas pixel⁻¹ in the distortion-free frame of Anderson & King 2006) and the time baseline between the two epochs (8.96 yr; Table 1). This yields transverse relative motions in mas yr⁻¹ in a frame closely aligned with the Galactic coordinate system.

Without absolute reference-frame tracers in this crowded field (e.g., Yelda et al. 2010; Sohn et al. 2012), we work exclusively with relative proper motions. Zero proper motion ${{\boldsymbol{\mu }}}_{f,0}$ is defined as the median observed rate of positional change for bulge objects across the entire field of view, without any selection for metallicity. The sample defining this proper-motion reference consists of stars that are not saturated in the deep exposures (these objects are at the bulge MS turnoff and fainter).¹³

Ca15 conducted extensive artificial-star tests to estimate measurement uncertainty in the proper motions, with artificial objects injected with proper motions into individual measurement frames to characterize the random proper-motion uncertainty as a function of apparent magnitude. Including random measurement uncertainty, random intrinsic uncertainty due to tracer star motion, and estimated residual distortion, the proper-motion uncertainty per coordinate is approximately ≲0.12 mas yr⁻¹ over the apparent magnitude range of interest (see Appendix A for details), easily sufficient to measure relative stellar motions in this field.

The result is a set of 339,193 objects with ACS/WFC positions, apparent magnitudes, and proper-motion estimates, all with uncertainties characterized as a function of apparent magnitude. Exploitation of these data is presented in Calamida et al. (2014, 2015) and Kains et al. (2017).

The WFC3 Bulge Treasury Project (BTS; program GO-11664; PI T. M. Brown) visited four fields in the bulge, with WFC3, including the SWEEPS field. The observations are described in detail in Brown et al. (2010); here we briefly summarize the characteristics relevant for the present paper.

In each field, observations were taken in UVIS/F390W (11,180 s), UVIS/F555W (2283 s), UVIS/F814W (2143 s), IR/F110W (1255 s), and IR/F160W (1638 s), with IR images (field of view 123'' × 136'') dithered in order to fully cover the UVIS observations (field of view 162'' × 162''). Good overlap was achieved with the SWEEPS ACS/WFC observations; nearly all the BTS objects in this field also fall within the SWEEPS ACS/WFC field of view (Figure 1).

Version 1 of the BTS catalog,¹⁴ which we use here, employed photometry and positions measured with daophotII (Stetson 1987; Brown et al. 2010). The resulting BTS v1 catalog lists 400,424 objects in the Sagittarius window with reported apparent magnitude in any of the BTS filters. Of these, 52,596 have measurements in all five of the BTS filters that are required to construct $[t],[m]$ estimates.

The BTS and SWEEPS catalogs were first cross-matched by equatorial coordinates. Although the absolute pointing of HST is accurate only to ∼01 (Gonzaga et al. 2012), with F814W observations in both data sets,¹⁵ matching of similar objects in both catalogs is straightforward (using F555W and F606W measurements in WFC3 and ACS/WFC, respectively, to refine the matches). For the first round of matching, a kd-tree approach was used to cross-match on the sphere, with a 5-pixel radius used for initial matching. In the second round, pixel positions in the two catalogs were cross-matched and fit using a general linear transformation for objects in the $18\leqslant {\rm{F}}814{\rm{W}}\leqslant 26$ range. While the population of good matches transitions to a background of mismatched objects at a radius of ∼2 pixels and larger, the vast majority of cross-matches were somewhat better, falling within a 1-pixel matching distance. The matching process resulted in a list of 47,537 objects with proper motions and seven-filter apparent magnitudes, with uncertainty estimates for all quantities.

Two aspects of the sample selection are worth highlighting. First, the selection region in the (F606W, F814W) CMD was chosen to be well clear of the MS turnoff, subgiant branch, and giant branch, to encompass as many stars as possible with good proper-motion measurements, and finally to capture a region over which the MS for a given population is reasonably free of curvature in the CMD. This selection region is shown in Table 2 and Figure 2. Second, the photometric metallicity and temperature indices include coefficients that amplify measurement uncertainty (particularly F110W and F160W, which appear in the temperature index $[t]$). For this reason, objects were only selected for further study for which all apparent magnitude uncertainties in the photometric catalog are smaller than 0.1 mag.

The successive selection steps isolating the sample for further study are detailed in Table 3. Of an initial sample of 339,193 SWEEPS objects and 400,424 BTS objects, 9700 (∼2.9%) were retained for further analysis.

Zoom In Zoom Out Reset image size

Table 3.  Selection Steps Used to Isolate the Proper-motion Sample for Further Study

Selection	N(remaining)	N(removed)
SWEEPS sample (Calamida et al. 2014)	339,193	⋯
Cross-matched with BTS	55,666	283,527
BTS measurements in all filters	47,537	8129
Within SWEEPS CMD selection region	10,225	37,312
σmag(ACS/WFC) < 0.1 mag	10,222	3
σmag(WFC3/UVIS) < 0.1 mag	10,209	13
σmag(WFC3/IR) < 0.1 mag	10,145	64
Clipping far outliers in $[t],[m]$	9,700	445

Note.  The cuts are cumulative, reading from top to bottom. The third step includes selecting out any rows for which any of the seven phometry and two proper-motion measurements are listed as a "bad" value in either the SWEEPS or BTS catalogs. In practice, this limits the sample to 18.5 ≤ F606W ≤ 27.5. The SWEEPS CMD selection region is shown in Figure 2. For the three instrumental configurations listed, objects must show photometric uncertainty <0.1 mag in all relevant filters: (F606W, F814W) for ACS/WFC, (F390W, F555W, F814W) for WFC3/UVIS, and (F110W, F160W) for WFC3/IR. Objects passing $[t],[m]$ clipping satisfy −3.50 ≤ $[t]$ ≤ −1.00 and −0.60 ≤ $[m]$ ≤ 0.40. See Section 3.2 for discussion.

Download table as:  ASCIITypeset image

The photometric indices $[t],[m]$ take the following form (Brown et al. 2009):

$$\begin{eqnarray}[t] & \equiv & (V-J)-\alpha (J-H)\\ [m] & \equiv & (C-V)-\beta (V-I),\end{eqnarray} \tag{ 1 }$$

with $\alpha \equiv E({\rm{F}}555{\rm{W}}-{\rm{F}}110{\rm{W}})/E({\rm{F}}110{\rm{W}}-{\rm{F}}160{\rm{W}})$ and $\beta \equiv E({\rm{F}}390{\rm{W}}\,-\,{\rm{F}}555{\rm{W}})/E({\rm{F}}555{\rm{W}}\,-\,{\rm{F}}814{\rm{W}})$, all of which have a dependence on stellar parameters. The median values of these stellar parameters for the proper-motion sample (${T}_{\mathrm{eff}}$ ≈ 4800 K and log(g) ≈ 4.6) were estimated from an isochrone chosen to overlap the observed sample (see Figure 2; several combinations of metallicity, age, and extinction were tried, indicating that the parameter range for this sample is roughly 4200 K ≲ ${T}_{\mathrm{eff}}$ ≲ 5200 K and $4.5\,\lesssim \mathrm{log}(g)\lesssim 4.7$).

The factors α, β are three-filter extinction ratios (Brown et al. 2009). Synthetic photometry was used to estimate the relationship between reddening and extinction for the objects of interest and to generate reddening vectors in the various filter combinations of interest. For a range of $E(B-V)$ values, pysynphot was used to generate synthetic stellar spectra, and the run of _AX against $E(B-V)$ was fit as ${A}_{X}={k}_{X}E(B-V)$ separately for all seven filters used in this study, over the range 0.0 ≤ $E(B-V)$ ≤ 1.5. The calculation was performed for ${T}_{\mathrm{eff}}$, $\mathrm{log}(g)$ appropriate to the SWEEPS CMD region chosen for proper-motion study (Figure 2). The process was repeated for low- and high-metallicity objects to estimate sensitivity of the extinction prescription to metallicity variation within the sample selected for further study, and for (${T}_{\mathrm{eff}}$, $\mathrm{log}(g)$) for objects at the median, minimum, and maximum ${T}_{\mathrm{eff}}$ within this sample to estimate spread of α, β along the sample.

This procedure requires a prescription for the extinction law toward the bulge. This extinction law appears to be somewhat nonstandard and strongly spatially variable, with some doubt in the literature about whether a single-parameter model can accurately reproduce observed behavior from the visible to the near-infrared (e.g., Nataf et al. 2016, and references therein). As the $[t],[m]$ indices use photometry over a very broad wavelength range (CVIJH, or λ ≈ 350–1700 nm), systematic uncertainties in the extinction prescription will in turn impact any inferences about the underlying metallicity distribution (this is one reason why we use $[t],[m]$ only to classify objects by relative [Fe/H] estimates).

To make progress, we adopted a single-parameter reddening law, but with ratio of selective to total extinction R_V = 2.5, as suggested by the investigations of Nataf et al. (2013).¹⁶ As this value is not among the standard parameterizations available in pysynphot, the coefficients ${A}_{X}/E(B-V)$ for the seven filters were estimated for R_V = 2.1 and R_V = 3.1 and linearly interpolated to R_V = 2.5.

Table 4 shows the k_X estimates for each filter, along with the coefficients α, β in the $[t],[m]$ indices. These are quite different from the MS coefficients reported in Brown et al. (2009), as expected since here we are targeting a specific population some way beneath the MS turnoff and have used a different prescription for extinction.

For a given choice of R_V the variation of all extinction-relevant quantities appears to be small within the sample of interest; α, β each vary by <0.1 between the two abundance sets tested and, for a given abundance, by ≲0.02 across the ${T}_{\mathrm{eff}}$ range of this sample. We adopt $(\alpha ,\beta )=(6.44,1.10)$ for the rest of this work.

The resulting ($[t],[m]$) distribution of objects is shown in Figure 3. Two concentrations are apparent: one near ($[t]$, $[m]$) = (−2.0, 0.15), and a second, more elongated concentration with major axis angled at about −45° in Figure 3, centered near ($[t]$, $[m]$) ≈ (−2.2, −0.1).

Zoom In Zoom Out Reset image size

To classify objects by relative metallicity and thus draw "metal-rich" and "metal-poor" samples for further study, the population highlighted in Figure 2 was characterized as a Gaussian Mixture Model (GMM) in ($[t]$, $[m]$) space, and members of the "metal-poor" and "metal-rich" samples were identified by their formal membership probability ${{\rm{W}}}_{{ik}}$ (see Appendix B). The number K of mixture components to use was determined by increasing K until the characterization stopped improving (see Appendix B.2 for details). At least two components seem to be required, but a four-component mixture model appears to provide the best representation of the $[t],[m]$ distribution.

We therefore adopt a four-component GMM to characterize the observed distribution in $[t],[m]$ space for the rest of this work. Table 5 shows the GMM parameters, while Figure 4 presents the model components visually. The two most significant components correspond roughly to visually apparent concentrations in Figure 3, together accounting for 91% of the mixture; these form our "metal-rich" and "metal-poor" samples. The remaining two components, making up about 6% and 3%, do not correspond to any physically obvious population. These two components might represent populations of outlier objects, or structure in the background in ($[t],[m]$). We retain these low-level components in the GMM for all subsequent work using the BTS catalog, but we do not interpret them as representing any intrinsic population component.

Zoom In Zoom Out Reset image size

Table 5.  Parameters of the Gaussian Mixture Model in $[t],[m]$ Space for Stars beneath the Main Sequence Selected for Further Study

k	Name	αk	${[t]}_{0}$	${[m]}_{0}$	${\sigma }_{[t][t]}^{2}$	${\sigma }_{[m][m]}^{2}$	${\sigma }_{[t][m]}^{2}$
			($\mathrm{mag}$)	(mag)	(mag2)	(mag2)	(mag2)
0	"Metal-poor"	0.557	−2.18	−0.09	0.0479	0.0143	−0.00742
1	"Metal-rich"	0.358	−1.97	0.15	0.0384	0.0043	−0.00187
2	⋯	0.026	−1.34	−0.05	0.0153	0.0397	0.00886
3	⋯	0.059	−2.87	0.01	0.0573	0.0255	0.00240

Note.  Reading left to right, columns indicate the component index k, its label (if any), its (rounded) mixture fraction α_k, the two components of its centroid, and the three unique components of the covariance matrix ${{\boldsymbol{V}}}_{k}$. See Section 3.5.

Download table as:  ASCIITypeset image

This four-component GMM provides the basis for our classification of objects by relative metallicity, with "metal-rich" and "metal-poor" objects corresponding to the two most significant components of the GMM (Table 5).¹⁷

A rough estimate for the centroid [Fe/H] values of the two samples may be drawn by charting [Fe/H] contours in the $[t],[m]$ diagram for synthetic stellar populations and interpolating to estimate [Fe/H] at the $[t],[m]$ locations of the corresponding GMM component centroids (see Appendix F.1 for more details on the synthetic stellar populations used). The GMM component centroids presented in Table 5 correspond to [Fe/H]₀ ≈ +0.18 for the "metal-rich" sample (using scaled-to-solar isochrones) and ${[\mathrm{Fe}/{\rm{H}}]}_{0}\approx -0.24$ for the "metal-poor" sample (using α-enhanced isochrones for this model component). These centroids are roughly consistent with values suggested from spectroscopic surveys (e.g., Hill et al. 2011; Zoccali et al. 2017).

For an object to be classified with the "metal-rich" or "metal-poor" sample, it must show formal membership probability ${{\rm{W}}}_{{ik}}\geqslant 0.8$ (see Equation (3); note that an object need not be classified with either sample when there are four model components). The shading in Figure 4 visualizes the membership probabilities ${{\rm{W}}}_{{ik}}$ associated with each mixture component. The threshold ${{\rm{W}}}_{{ik}}\geqslant 0.8$ was chosen as a trade-off between sample purity (typical objects should not fall into more than one model component at the chosen threshold) and the need to have a sufficient sample size (at least a few thousand) to permit the dissection of the proper motions by relative photometric parallax with sufficient resolution to chart the rotation curves.

Assigning relative photometric parallax ($\pi ^{\prime}$) is the final step required before proper-motion rotation curves can be charted, with reference to fiducial ridgelines for the "metal-rich" and "metal-poor" samples. The fiducial ridgelines themselves were determined by a simple empirical fit to the density of each sample in the SWEEPS CMD. A second-order polynomial adequately represents the median samples and allows very rapid evaluation of relative photometric parallax. Figure 5 shows the adopted fiducial ridgelines for the "metal-rich" and "metal-poor" samples  in the SWEEPS CMD, while their parameters are given in Table 6.

Zoom In Zoom Out Reset image size

Table 6.  Ridgeline Parameters in the SWEEPS Color–Magnitude Diagram, for the "Metal-poor" and "Metal-rich" Samples

k	Name	a0	a1	a2
		($\mathrm{mag}$)		(mag−1)
0	"Metal-poor"	−16.906	48.557	−14.954
1	"Metal-rich"	−5.720	33.458	−10.293

Note.  These purely empirical ridgelines are used to rapidly evaluate photometric parallax for objects in each sample and take the form ${\rm{F}}814{\rm{W}}={{\rm{\Sigma }}}_{j}{a}_{j}{x}^{j}$, with x the (F606W – F814W) color. See Section 3.5 for discussion.

Download table as:  ASCIITypeset image

Having drawn "metal-rich" and "metal-poor" samples from the BTS photometry, along with fiducial sequences in the SWEEPS CMD for the two samples, the "metal-rich" and "metal-poor" rotation curves can be charted. Figure 6 shows the raw distribution of longitudinal proper motion μ_l and relative photometric parallax for the "metal-rich" and "metal-poor" samples, with trends presented in Figure 7.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

All proper motions in this work were measured relative to the same proper-motion zero-point, defined without reference to any selection by metallicity (Section 2.1). To the extent that the "metal-rich" and "metal-poor" samples trace bulge objects with different spatial distribution and/or kinematic behavior, however, the average proper motions of bulge objects in the two samples might differ.

We therefore estimated the proper motion corresponding to the fiducial sequence for each sample. For this "central" proper motion, we used the median proper motion over those sample members with relative photometric parallax between ${\rm{\Delta }}{\pi }_{-}^{{\prime} }$ and ${\rm{\Delta }}{\pi }_{+}^{{\prime} }$ magnitudes nearer to and farther than the fiducial, respectively.¹⁸ We adopted ${\rm{\Delta }}{\pi }_{+}^{{\prime} }$ = 0.05, corresponding roughly to Δ$D$ ≈ 0.18 kpc at the distance of the bulge (Section 4.1). The central proper motion for the "metal-rich" sample is then ${({\mu }_{l},{\mu }_{b})}_{\mathrm{MR}}^{0}=(+0.019,+0.19)$ mas yr⁻¹ from 381 surviving objects, while for the "metal-poor" sample we found ${({\mu }_{l},{\mu }_{b})}_{\mathrm{MP}}^{0}=(-0.12,+0.32)$ mas yr⁻¹ from 290 surviving objects. The proper motion corresponding to the "metal-rich" fiducial was thus found to be offset from that of the "metal-poor" fiducial by about ${({\mu }_{l},{\mu }_{b})}_{\mathrm{MR}\mbox{--}\mathrm{MP}}^{0}$ ≈ $(+0.14,\,-0.13)$ mas yr⁻¹.

With a difference in rotation curves suggested from the behavior of μ_l against relative photometric parallax, the next step is to chart the distance variation of the (l, b) proper-motion ellipse. The approach shares several similarities to that reported in Cl08; relative photometric parallaxes were assigned to each star with reference to the fiducial sequence (appropriate for the metallicity sample with which the star was identified) and the sample partitioned into bins of relative photometric parallax $\pi ^{\prime}$, with bin widths adjusted so that each bin has the same number of objects.

The proper-motion distribution within each bin was fit as a 2D Gaussian, with centroid proper motion ${{\boldsymbol{\mu }}}_{0}$ and covariance matrix ${{\boldsymbol{V}}}_{\mu }$. Uncertainties in fitted quantities were estimated by parametric bootstrapping: synthetic samples for each bin were drawn from the best-fit model, perturbed by the estimated proper-motion uncertainty, and the distribution of recovered parameters over the bootstrap trials was adopted as the estimated parameter uncertainties. Because this process can be sensitive to outliers, a single pass of sigma-clipping was applied to the proper-motion sample within each distance bin using a ±3σ threshold; this typically removed roughly 1%–2% of the points per bin, with the exeption of the most distant $\pi ^{\prime}$ bin (see Tables 16 and 17 in Appendix I).

Several improvements were made over the analysis reported in Cl08. For example, rather than subtracting the estimated proper-motion uncertainty in quadrature from the model covariances after fitting, the "extreme-deconvolution" formulation of Bovy et al. (2011) was used, which incorporates estimated measurement uncertainty as part of the fitting process (see Appendix B). We experimented with a multicomponent GMM within each $\pi ^{\prime}$ bin for each sample but found a single component adequate (see also Section 5.8). The estimates of proper-motion uncertainty themselves have also been improved compared to Cl08, both in the characterization of random uncertainty through the artificial-star tests of Ca15 and through improved characterization of residual relative distortion (Kains et al. 2017). Details of the adopted uncertainty estimates are presented in Appendix A; for the apparent magnitude range of interest, the total proper-motion uncertainty estimates (${\epsilon }_{i}\lesssim 0.12$ mas yr⁻¹) are much smaller than the intrinsic proper-motion dispersion of the bulge (∼3 mas yr⁻¹).

Figure 11 presents the rotation and dispersion curves of the "metal-rich" and "metal-poor" samples expressed in terms of $(D,v)$. The conversion of these quantities from the measured ($\pi ^{\prime} ,\mu$) requires the reference distance ${D}_{0}$ corresponding to the fiducial sequences for the two samples. The reference distance was set by taking literally the distance modulus (m − M)₀ = 14.45 suggested by studies of the SWEEPS CMD (Ca14), which in turn suggests reference distance (${D}_{0}=7.76\,\mathrm{kpc}$). We assigned this reference distance to both the "metal-rich" and "metal-poor" samples (a choice we examine critically in Section 5.4).

Consistent with the simple treatment in Figure 7 and Section 3.6, the "metal-rich" sample shows a higher-amplitude rotation curve than does the "metal-poor" sample, with both a steeper slope and about a factor of ∼2 greater difference in mean transverse velocity $\langle {v}_{l}\rangle$ between the near side and far side of the bulge than for the "metal-poor" sample.

To quantify the rotation curve discrepancies between the samples, a simple straight-line model was fit to their rotation curves for distances close to the fiducial for each sequence. This interval was estimated separately for the two samples since their rotation curves appear to level off at different distances from the fiducial (Figure 11). For the "metal-rich" sample the gradient was estimated over the interval ${D}_{0}\pm 0.80\,\mathrm{kpc}$ in the ($D$, ${v}_{l}$) curve (corresponding to −0.24 mag ≲ $\pi ^{\prime}$ ≲ +0.21 mag). The rotation curve of the "metal-poor" sample remains sloped over a broader range, so the fitting interval ${D}_{0}$ ± 1.4 kpc was used (so −0.44 mag ≲ $\pi ^{\prime}$ ≲ +0.36 mag). For both samples the rotation curve amplitude was estimated from the intervals where the rotation curves level off, covering two to three bins each outside the sloped region (Figure 12 indicates the regions used to estimate the rotation curve slopes and amplitudes). The 1σ ranges of $\langle {\mu }_{l}\rangle$ and $\langle {v}_{l}\rangle$ from the parametric bootstrap trials were used as estimates of measurement uncertainty in each distance bin, and the trends were fitted to each of the ($\pi ^{\prime}$, $\langle {\mu }_{l}\rangle$) and ($D$, $\langle {v}_{l}\rangle$) rotation curves separately (rather than transforming the proper-motion trends into velocity trends after fitting). We did not attempt to deproject velocities to circular speeds (as discussed in Cl08) but merely attempted to characterize observed trends.

Figure 12 and Table 7 show the results. The ratio of the gradients $B$ was found to be ${({B}_{\mathrm{MR}}/{B}_{\mathrm{MP}})}_{l}$ = $3.70\pm 0.68$, while the ratio of amplitudes $A$ is ${({A}_{\mathrm{MR}}/{A}_{\mathrm{MP}})}_{l}$ = $2.29\pm 0.35$. Thus, a ratio in rotation curve slopes was detected at approximately $5.4\sigma$, while for the velocity amplitude the ratio was detected at roughly $6.5\sigma$.

The rotation curves in Galactic latitude (Figure 8, bottom panel) visually suggest gentle trends from near side to far side, consistent with previous measurements (e.g., Cl08, Soto et al. 2014). Table 8 reports the straight-line characterization of the Galactic latitude rotation curves, where here we characterized the trend as a straight-line fit within ±2.0 kpc from the fiducial distance ${D}_{0}$. The behaviors of the two samples in μ_b are statistically similar, with gradient ratio ${({B}_{\mathrm{MR}}/{B}_{\mathrm{MP}})}_{b}$ = $0.90\pm 0.87$. We do not consider this to represent a secure detection of differing rotation curves in the direction of Galactic longitude.

The velocity dispersion profiles (measured as major- and minor-axis lengths of the proper motion and velocity ellipsoids; Figures 9 and 11) also show differences between the samples. Both samples show a broadly centrally peaked velocity dispersion pattern against line-of-sight distance (Figure 11), with the "metal-rich" sample showing a narrower peak, particularly in the major-axis dispersion. The "metal-rich" sample also shows generally lower velocity dispersion by ∼10%, particularly in terms of the velocity ellipse minor axis.

Zoom In Zoom Out Reset image size

Consistent with previous studies (e.g., Soto et al. 2014), the proper-motion ellipse appears to be weakly elongated, with the "metal-rich" population possibly the more elongated of the two samples (with axis ratio $b/a\approx 1.29\pm 0.05$ at $\pi ^{\prime}$ = 0 compared to $b/a\approx 1.13\pm 0.05;$ see the top panel of Figure 10). However, the two axis-ratio trends show considerable bin-to-bin scatter.

Zoom In Zoom Out Reset image size

The proper-motion ellipse major-axis position angle also shows trends with relative photometric parallax, although possibly at lower statistical significance than the trends reported in Cl08 despite a much longer time baseline for proper motions (Ca14). This reduced significance may be due to the reduced sample size admitted by the cuts in $[t],[m]$ employed in this work. It may be that only the "metal-rich" sample substantially shows the proper-motion ellipse tilt with distance, with position angle rising to the 20°–40° range (this tilt is strongly influenced by projection effects; see Section 5.1 and particularly Equation (2) of Cl08). Because the "metal-poor" population tends to be less elongated, its position angle trends are also detected at lower significance.

The very nearest relative photometric parallax bins show behavior consistent with a foreground population dominated by Galactic rotation. This seems particularly clear for the "metal-rich" sample, which shows a much more strongly elongated proper-motion ellipse for the nearest bin ($a/b\approx 2.0\pm 0.11$) and position angle consistent with zero (consistent with differential rotation in Galactic latitude).

Differences in apparent magnitude distribution other than those due to distance spread would contribute to differences in the inferred $\pi ^{\prime}$ distributions for the "metal-rich" and "metal-poor" samples. If sufficiently severe, this differential blurring in $\pi ^{\prime}$ might cause two intrinsically identical rotation curves to be erroneously measured as discrepant. In the sense of our findings, the "metal-poor" sample might be artificially blurred in $\pi ^{\prime}$ compared to the "metal-rich" sample, which would produce an apparent rotation curve discrepancy where none were present.

Several phenomena might lead the "metal-poor" sample to exhibit greater apparent magnitude scatter than the "metal-rich" sample. First, since the "metal-poor" ridgeline in the SWEEPS CMD is slightly fainter than the "metal-rich" ridgeline, the "metal-poor" objects may be subject to increased photometric uncertainty. Second, at least in principle, if the extinction distribution experienced by the two samples were in some way different, this could lead to a broader apparent magnitude distribution for the "metal-poor" sample. Third, differences in binary fraction between the samples might cause the relative photometric parallax distribution of the two samples to differ, although the nature, magnitude, and direction of such effects may be complex and indeed depend on the class of binaries probed (e.g., Gao et al. 2014).

Finally, differences in the intrinsic photometric scatter between the "metal-rich" and "metal-poor" samples might amplify differences between the rotation curves. Our own VLT spectroscopy, as well as spectroscopic campaigns from the literature (e.g., Hill et al. 2011; Zoccali et al. 2017), suggests that the [Fe/H] spread for the "metal-poor" population is greater than for the "metal-rich" population, which would in turn contribute greater $\pi ^{\prime}$ scatter in the "metal-poor" population.

We have performed simple Monte Carlo tests to determine whether perturbations in the inferred distance distribution can be responsible for the differences in rotation curves between "metal-rich" and "metal-poor" samples. Appendices D and E provide details.

In the course of investigating the impact of the differential [Fe/H] distribution on the $\pi ^{\prime}$ distribution, it became apparent that the BaSTI set of artificial stellar population methods used to generate synthetic [Fe/H] distributions was (at the time of this work) imposing an apparently artificial population truncation. Appendix F provides details, with the method we adopted to mitigate this selection effect discussed in Appendix E.

Perturbations were tested as a result of additional photometric uncertainty or differential extinction variations (Appendix D.1), differences in the fraction of unresolved binaries (Appendix D.2), and differences in the photometric parallax spread caused by differing intrinsic spreads in metallicity (Appendix E). In all scenarios, either the effect is too small to bring the rotation curves into agreement (for binaries), or the required perturbation is too large to have gone unnoticed in previous studies (for extinction), possibly by an order of magnitude (for photometric uncertainty). The strongest single contributor of relative photometric parallax mixing is intrinsic difference in metallicity spread between the samples; this likely contributes differential distance mixing up to one-third the amount required to artificially reproduce the observed discrepancy in rotation curves. Since independent sources of additional photometric scatter would presumably add in quadrature, their combination is very unlikely to be sufficient to bring about the observed discrepancies in trends.

We therefore conclude that differential distance scatter is not responsible for the difference in rotation curves or $\pi ^{\prime}$ distributions, due to additional photometric uncertainty, differential extinction, differences in the unresolved binary populations, or differences in metallicity spread between samples.

While blurring in relative photometric parallax would tend to artificially increase the difference between trends in the "metal-rich" and "metal-poor" samples, cross-contamination of the samples in ($[t],[m]$) would tend to artificially reduce these differences. While we have used reasonably conservative thresholds in drawing our "metal-rich" and "metal-poor" samples, genuinely metal-rich objects might be moved into the "metal-poor" sample by measurement uncertainty, and vice versa.

Because of the complexities involved in rigorous reconstruction of the observed distributions (e.g., Gennaro et al. 2015), full exploration of this cross-contamination is deferred to future work. We have performed a simple Monte Carlo contamination test for the formal membership probability threshold ${{\rm{W}}}_{{ik}}$ > $0.8$ used in this work (Appendix G). Under the assumptions of that test, we find that the "metal-rich" sample is contaminated at the ∼5% level (mostly from the "metal-poor" sample), while the "metal-poor" sample is contaminated at the ≲1% level (mostly due to the "metal-rich" sample, but with some contribution from background component k = 3 in Table 5). This is not severe enough for the observed low-amplitude "metal-poor" rotation curve to be due to sample contamination from a small population of objects following the kinematics of the "metal-rich" sample.

Our treatment of the impact of extinction on the photometry (and therefore the relative photometric parallax) assumes that the extinction is constant over the line-of-sight distances of interest. Violations of this assumption might in principle influence the trends we observe, by artificially broadening the line-of-sight distribution (with stars more affected by extinction appearing farther from the observer).¹⁹ Here we examine the likely impact on our main results of extinction variations along the line of sight.

A few studies have mapped the 3D extinction distribution out to the far side of the bulge (e.g., Marshall et al. 2006; Schultheis et al. 2014). Particularly for sight lines close to the Sagittarius window we study here, most of the extinction at these distances takes place at distances $D$ ≲ 5 kpc from the Sun, possibly broken into two foreground concentrations (at $D$ ≈ 3 and 5 kpc; e.g., Marshall et al. 2006). Thus, along our sight line, extinction variations within the bulge are likely to be small compared to variations in the foreground disk, an interpretation consistent with the photometry of red clump stars in this field (Cl08). Hence, the trends we find for line-of-sight distances 5 kpc ≲ $D$ ≲ 11 kpc—the main sample of interest—are likely unaffected.

By symmetry we might expect additional extinction from at least one dust screen at distances $D$ ≳ 11 kpc owing to spiral structure on the far side of the bulge (see, e.g., Figure 10 of Schultheis et al. 2014). This would be mitigated somewhat by the slightly tilted path of our line of sight compared to the Galactic plane; at Galactic latitude b_J2000.0 = −265, our sight line is already ≈480 pc below the Galactic midplane when it reaches $D$ = 11 kpc, roughly where it might intersect the first dust concentration on the far side of the bulge (compared to ≈210 pc at 5 kpc), suggesting that far-side extinction may likely be somewhat weaker than experienced in the foreground.

We therefore conclude that indeed the photometric parallaxes for objects closer than $D$ ≈ 5 kpc and farther than $D\approx 11\,\mathrm{kpc}$ may have been assigned photometric parallaxes that are artificially close and far, respectively. However, as those distances are outside our main regions of interest, this does not impact any of the trends that we report.

All proper motions in this work are reported relative to the same proper-motion zero-point ${{\boldsymbol{\mu }}}_{{\rm{F}},0}$, which is defined as the average proper motion of astrometrically well-measured stars, whatever their metallicity (Section 2.1). We have chosen not to apply separate proper-motion zero-points to the two samples, for example, by forcing the two samples to each show ${\mu }_{l}$ = 0.0 mas yr⁻¹ at $\pi ^{\prime}$ = 0.0 mag, but opt to keep the proper motions in the same reference frame for each sample to allow direct comparison between ($\pi ^{\prime}$, μ) rotation curves.

We then find that the central proper motions for the "metal-rich" and "metal-poor" samples (i.e., the median proper motions for stars near their fiducial sequences) differ by ${({\mu }_{l},{\mu }_{b})}_{\mathrm{MR}-\mathrm{MP}}^{0}\approx (+0.14,\,-0.13)$ mas yr⁻¹ (Section 3.6); equivalently, the rotation curves do not meet at ($\pi ^{\prime}$, ${\mu }_{l}$) = (0 mag, 0 mas yr⁻¹).

These discrepancies could be due to differences in the mean intrinsic velocities of the fiducial stars between the samples, or differences in the line-of-sight distance at which the fiducual stars are found, or a combination of the two. Since the fiducial sequences for the two samples are determined from their CMD population densities (see Figure 5), their central proper motions could well differ if their densest observed regions occurred at different distances. This could occur naturally if the two samples are oriented differently in the Galactic plane (in which case the relationship between the mean velocity ${\boldsymbol{v}}$ and its transverse velocity components ${v}_{l}$ and ${v}_{b}$ would also differ between the samples).

If the fiducial stars for the two samples do indeed lie at different distances ${D}_{0,k}={D}_{0}+{{\rm{\Delta }}}_{k}$ (with Δ_k giving the distance offset for a particular sample), then the appropriate conversion from ($\pi ^{\prime}$, ${\mu }_{l}$) to ($D$, ${v}_{l}$) will also differ, in turn impacting the difference in rotation curve velocity amplitude for the two samples. Figure 13 shows how the velocity rotation curve is impacted by shifting the "metal-rich" fiducial by distance offset Δ kpc from the "metal-poor" fiducial. Bringing the "metal-rich" fiducial closer than the "metal-poor" one does reduce the velocity amplitude discrepancy between the two samples. However, to bring the velocity amplitudes of the two samples into agreement, the "metal-rich" fiducial would need to be brought closer by $| {\rm{\Delta }}| \gtrsim 2\,\mathrm{kpc}$, which seems unlikely for samples so close to the Galactic rotation axis (at $l\approx +1\buildrel{\circ}\over{.} 26$). Unless the spatial distributions of the two samples really are radically different, then, we consider it unlikely that a difference in fiducial distance between the samples can by itself produce the observed difference in velocity rotation curve amplitude we are measuring.

Zoom In Zoom Out Reset image size

While an offset Δ in fiducial distance scales the velocity amplitude by a corresponding amount, an offset ${\rm{\Delta }}{{\boldsymbol{\mu }}}_{{\rm{F}},0}$ in the proper-motion zero-point produces a systematic shift ${\rm{\Delta }}(d{\boldsymbol{v}}/{dD})=-4.74{\rm{\Delta }}{{\boldsymbol{\mu }}}_{{\rm{F}},0}$ in the velocity gradient. If ${{\boldsymbol{\mu }}}_{{\rm{F}},0}$ (the average proper motion of well-measured bulge stars of all metallicities) and ${D}_{0}$ (the average distance to bulge stars of all metallicities) are both determined from the same set of stars, then we would have ${\rm{\Delta }}{{\boldsymbol{\mu }}}_{{\rm{F}},0}$ = 0 mas yr⁻¹ and thus ${\rm{\Delta }}(d{\boldsymbol{v}}/{dD})$ = 0 $\mathrm{km}\,{{\rm{s}}}^{-1}\,{\mathrm{kpc}}^{-1}$.

However, in reality the sets of stars used to estimate ${{\boldsymbol{\mu }}}_{{\rm{F}},0}$ and ${D}_{0}$ will in general differ. The fiducial distance ${D}_{0}$ is estimated from the distribution of "extreme-bulge"(EB) stars showing ${\mu }_{l}\lt -2.0$ mas yr⁻¹ (e.g., Ca15), while ${{\boldsymbol{\mu }}}_{{\rm{F}},0}$ is estimated from stars below the MS turnoff without any proper-motion selection. Thus, although stars well measured astrometrically tend also to be well measured photometrically, differing selection effects in the determination of ${D}_{0}$ and ${{\boldsymbol{\mu }}}_{{\rm{F}},0}$ will still lead to a global systematic offset in ${\rm{\Delta }}{{\boldsymbol{\mu }}}_{{\rm{F}},0}$.

The true value of ${\rm{\Delta }}{{\boldsymbol{\mu }}}_{{\rm{F}},0}$ is unknown; however, we can form a rough estimate as the median proper motion of the population traced by the EB objects used to estimate ${D}_{0}$ (Ca14). The EB objects have highly negative proper motions by construction, but we can estimate the median proper motion of the underlying population that they trace by estimating the median photometric parallax $\langle \pi ^{\prime} \rangle$ of the EB tracers and applying the proper-motion rotation curve characterization ${\mu }_{l}(\pi ^{\prime} )$ of Section 4.1. Since we have performed this characterization separately for "metal-rich" and "metal-poor" samples, we can estimate ${\rm{\Delta }}{{\boldsymbol{\mu }}}_{{\rm{F}},0}$ separately from the two samples. Applying the kinematic cut ${\mu }_{l}$ ≤ −2.0 mas yr⁻¹ to extract EB tracers for the "metal-rich" and "metal-poor" populations, we find median photometric parallax $\langle \pi ^{\prime} {\rangle }_{\mathrm{EB},\mathrm{MR}}$ ≈ +0.025 mag and $\langle \pi ^{\prime} {\rangle }_{\mathrm{EB},\mathrm{MP}}$ ≈ −0.002 mag for EB objects in the "metal-rich" and "metal-poor" samples, respectively (so that objects with ${\mu }_{l}$ = 0 mas yr⁻¹ lie slightly in front of the EB population, as expected; see Section 5.9). This suggests that the underlying population traced by the EB objects—corresponding to the fiducial distance ${D}_{0}$—has median longitudinal proper motion $\langle {\mu }_{l}{\rangle }_{\mathrm{MR}}\approx -0.17\,\mathrm{and}\,\langle {\mu }_{l}{\rangle }_{\mathrm{MP}}\approx 0.00$ mas yr⁻¹. These figures are likely also sensitive to differing intrinsic proper motion and distance distributions between the samples, but this estimate suggests that the proper-motion difference between the sample from which the proper-motion reference frame was set and the sample from which ${D}_{0}$ was estimated is not larger than ${\rm{\Delta }}{{\boldsymbol{\mu }}}_{{\rm{F}},0}$ ≲ 0.17 mas yr⁻¹. Thus, the systematic velocity gradient uncertainty ${\rm{\Delta }}(d{\boldsymbol{v}}/{dD})$ may be on the order of ∼$0.8\,\mathrm{km}\,{{\rm{s}}}^{-1}\,{\mathrm{kpc}}^{-1}$.

Systematic uncertainty in the proper-motion zero-point ${{\boldsymbol{\mu }}}_{{\rm{F}},0}$ may therefore impact the ratio of longitudinal velocity gradients reported in Section 4.1 by ≲10% (Table 7), which is too small to materially affect the main results or conclusions we report.

We are finally in a position to answer the question posed by Section 1.1. Our "metal-rich" and "metal-poor" rotation curves are inconsistent with each other at $\sim 5.4\sigma$ for the rotation curve slope and $\approx 6.5\sigma$ for the near-side−far-side rotation amplitude (Section 4.1).

These proper-motion-based results stand in strong contrast to determinations of the radial velocity rotation trends, which either find weak if any discrepancy between "metal-rich" and "metal-poor" samples (e.g., Ness et al. 2013b; Williams et al. 2016; Zoccali et al. 2017) or require a large contribution from samples with [Fe/H] < −1.0 to produce a discrepancy in rotation curves (e.g., Minniti 1996; Kunder et al. 2016; note that the fraction of stars in our sight line with [Fe/H] < −1.0 is likely low; Zoccali et al. 2017). It seems unlikely that this discrepancy between our proper-motion-based study and these radial-velocity-based studies can be due purely to any differences between the uses of giant and dwarf stars as tracers, since microlensed dwarf stars also show no strong differences in mean radial velocity between metal-poor and metal-rich stars (or, for that matter, between stars younger and older than ∼7 Gyr; Bensby et al. 2017).²⁰

The most likely explanation for the difference between proper-motion and radial velocity results is the strong difference in the way the studies sample the inner Milky Way. Detailed comparison of our new observational indications with model prediction is deferred to future work; however, a likely scenario to explain the differences can be outlined as follows. In the simulations of Debattista et al. (2017), the final orbital configuration of bulge stars depends on their (radial) velocity dispersion before bar formation. To the extent that [Fe/H] and radial velocity dispersion correlate with each other (or, equivalently, each correlate with time of formation), the "kinematic fractionation" resulting from bar formation might well leave metal-rich stars with a higher fraction of elongated orbits than for metal-poor objects. That the "X" structure is observed to preferentially contain metal-rich objects (e.g., Vásquez et al. 2013) supports the notion that stars within our "metal-rich" and "metal-poor" samples on average move along differently shaped orbits, while the spatial structures traced by Mira variables of different pulsation period ranges suggest that samples with differing spatial configuration can coexist in the same volume at the present day (Catchpole et al. 2016). The "metal-rich" and "metal-poor" samples may then show quite different transverse velocity distributions as a function of line-of-sight distance, even if the distributions produce similar mean velocities when averaged along the line of sight owing to radial velocity survey selection effects. In this scenario, only by dissecting the population by line-of-sight distance (or its proxy, $\pi ^{\prime}$) can the differing velocity distributions of the two coexisting samples be distinguished.

Since stars in the "metal-rich" and "metal-poor" samples all move through the same present-day potential, by detecting differences in the proper-motion-generated rotation curve, we may well be detecting differences in orbital anisotropies between metal-rich and metal-poor bulge objects. While detailed prediction is a topic of ongoing work, the differences we detect seem qualitatively reasonable at present.

Having shown that the proper-motion-based rotation curve does show discrepancy between "metal-rich" and "metal-poor" populations, the necessary next step is to extend our approach to more sight lines within the inner bulge. By comparing metallicity-dissected proper-motion-based rotation curves between fields, the trends with location in the bulge can be charted empirically, allowing a sharper test of the true variation of bulge rotation with the metallicity of the sample probed. This work is deferred to a future communication.

Both the "metal-rich" and "metal-poor" samples show a clear central peak in velocity ellipse major-axis length near the distance interval where the samples are the most densely populated (Figure 11). The peak persists in the velocity ellipse minor axis for the "metal-rich" sample but is rather less clear in the "metal-poor" sample. This is broadly similar to the trends found from the combined population in previous studies (e.g., Cl08, Soto et al. 2014). We note a rough, qualitative similarity with the curve of radial velocity against Galactic longitude for Galactic latitude b = −2° (see the middle left panel of Figure 12 of Zoccali et al. 2017); however, we remind the reader that the radial velocity and proper-motion trends cannot directly be compared because they suffer from differing selection effects. That the proper-motion dispersion of the "metal-poor" component is generally slightly larger than that of the "metal-rich" one (particularly along the minor axis) is qualitatively consistent with expectations that a metal-poor, less rotationally supported population should show higher velocity dispersion (e.g., Ness et al. 2013b; Debattista et al. 2017).

We may also be detecting the velocity dispersion "inversion" detected at the innermost fields in radial velocity studies (Babusiaux et al. 2014; Zoccali et al. 2017). Consistent with the low-latitude radial velocity dispersion trends, the proper-motion-based velocity dispersion might also be greater for the "metal-rich" sample than for the "metal-poor" sample at the distance bins closest to the center of the bulge (see Figure 11). For the innermost bulge regions, the proper-motion-based "metal-rich" velocity dispersions also show steeper gradient than the "metal-poor" ones, but with the gradient against line-of-sight distance rather than Galactic longitude, with the innermost distance bin possibly showing slightly greater velocity dispersion for the "metal-rich" sample.

The tendency of the "metal-poor" sample to show greater dispersion in relative photometric parallax (or, correspondingly, in distance $D$) is qualitatively consistent with the "kinematic fractionation" of Debattista et al. (2017). Under that mechanism, more "metal-poor" populations also initially had greater radial velocity dispersion, leading to a distinct (and broader) present-day spatial distribution when compared to the most "metal-rich" objects.

The difference we find in line-of-sight distribution between the "metal-rich" and "metal-poor" samples might also be consistent with the observations of Catchpole et al. (2016), who find differing bar angles and degrees of central concentration for Mira variables of different ages. However, the interplay between age and metallicity of bulge stars is likely not simple. For example, while a gentle relationship may exist between [Fe/H] and the fraction of stars younger than about 8 Gyr (e.g., Figure 14 of Bensby et al. 2017 or Figure 10 of Bernard et al. 2018), the microlensing spectroscopic surveys suggest that stars can take any metallicity value (for [Fe/H] ≳ −1.0) for any age. How the predictions of Catchpole et al. (2016) translate into predictions for the two samples here is deferred to future work.

In addition to any continuous metallicity–velocity correlation, the samples may include populations from distinct entities within the bulge region, whether interlopers from the halo (e.g., Koch et al. 2016), any thick-disk component (e.g., Ness et al. 2013a), or a small "classical" bulge component (Kunder et al. 2016).

We have therefore performed the exercise of decomposing each of the "metal-rich" and "metal-poor" samples into two-component GMMs (in ${\mu }_{l}$ only), to determine whether any minority component is distinguishable within the rotation curves formed from the two samples (Figures 14 and 15). No minority population is detected in either sample; indeed, when a two-component GMM is used, the two centroids track the mean rotation curve within each sample roughly symmetrically about the mean rotation curve, while each subcomponent has roughly equal weight in the mixture.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We therefore conclude that a minority component with discrepant rotation curve is not required in either the "metal-rich" or "metal-poor" sample, but due to the small sample size (≈2000 stars in total per sample), we cannot at this stage rule out its presence. Direct comparison with population models may allow upper limits to be set on the presence of any minority component within each sample, but this is deferred to future work.

Photometric studies of the bulge typically impose a condition ${\mu }_{l}\lt -2.0$ mas yr⁻¹ to isolate a clean-bulge sample for further study (e.g., Kuijken & Rich 2002; Calamida et al. 2014), although there are exceptions (e.g., Bernard et al. 2018).²¹ This procedure is appropriate because in the sight lines typically studied near the Galactic center the foreground disk population typically shows proper motion relative to the mean-bulge population of Δμ_l ≈ +4 mas yr⁻¹, as suggested by direct comparison of the proper motions of bulge giant branch stars with those of the upper MS population of (mostly) disk foreground stars (e.g., Cl08; Soto et al. 2014).

To investigate whether and how a simple cut on ${\mu }_{l}$ imposes selection effects on the two samples, we computed the sample counts, fractions, and volume densities for objects that would pass the longitudinal proper-motion cut (μ_l < −2.0 mas yr⁻¹). The results are plotted in Figure 16 and presented in tabular form in Tables 16 and 17 in Appendix I.

Zoom In Zoom Out Reset image size

We find that the cut (μ < −2.0 mas yr⁻¹) admits very few foreground objects from either sample; both samples show fewer than two objects passing this cut for the closest distances ($D\lesssim 4.53\,\mathrm{kpc}$ and ≲4.10 kpc for "metal-rich" and "metal-poor" samples, respectively; see Tables 16 and 17 in Appendix I), while the mean proper motion $\langle {\mu }_{l}\rangle$ of the foreground population climbs strongly for the closest distance bins (Figure 8). We therefore confirm that the traditional proper-motion cut (μ < −2.0 mas yr⁻¹) does indeed remove nearby objects cleanly for the SWEEPS field.²²

Beyond this, however, the dissection by relative abundance has revealed several interesting selection effects among the kinematically cleaned sample (Figure 16).

First, as expected, there is a bias toward the far side of the bulge, but this bias is much stronger in the "metal-rich" sample than for the "metal-poor" one; indeed, the fraction of "metal-poor" objects passing the kinematic cut is almost flat with inferred distance between ~5 and 9 kpc.

Second, the raw counts of sources thus isolated in the "metal-rich" and "metal-poor" samples are of similar orders of magnitude. Considering sample sizes that pass the kinematic cut at inferred distances between 6.4 and 9.1 kpc (chosen to encompass the bulge populations; see Tables 16 and 17), the total counts in each sample are 521 ± 19 and 507 ± 19 for the "metal-rich" and "metal-poor" samples, respectively (the uncertainties, estimated from the quadrature sum of parametric bootstrap uncertainty estimates in these counts for each bin, are almost certainly underestimates). With total sample sizes within this distance range of 2181 (1783) for the "metal-rich" ("metal-poor") samples, this translates into fractions 24% ± 1% (28% ± 1%) of the "metal-rich" ("metal-poor") samples that pass the kinematic cut. Thus, of objects in this distance range, the kinematic cut appears to slightly favor the "metal-poor" sample, although the difference is small.

In principle, a population of compact objects among the foreground population might fall into the farther distance bins for the "metal-poor" sample,²³ polluting a sample with bulge-like motions with a small population showing disk-like motions. However, with the foreground disk population at ∼10% of the total (Ca14) and with a substantial white dwarf population perhaps unlikely for a typical "young" foreground population, we do not consider this a significant contaminant, and we leave exploration of the impact of foreground white dwarfs to future work.

From the definition of the $[t],[m]$ indices (Equation (1)), uncertainty propagation produces an approximation for the appropriate measurement uncertainty covariance ${{\boldsymbol{S}}}_{i}$ for each data point, which we reproduce here for convenience. We adopt

$$\begin{eqnarray}{{\boldsymbol{S}}}_{i}={\left(\begin{array}{cc}{\sigma }_{t}^{2} & {\sigma }_{{tm}}^{2}\\ {\sigma }_{{mt}}^{2} & {\sigma }_{m}^{2}\end{array}\right)}_{i}\\ ={\left(\begin{array}{cc}{\sigma }_{V}^{2}+{(1+\alpha )}^{2}{\sigma }_{J}^{2}+{\alpha }^{2}{\sigma }_{H}^{2} & -(1+\beta ){\sigma }_{V}^{2}\\ -(1+\beta ){\sigma }_{V}^{2} & {(1+\beta )}^{2}{\sigma }_{V}^{2}+{\sigma }_{C}^{2}+{\beta }^{2}{\sigma }_{I}^{2}\end{array}\right)}_{i,},\end{eqnarray} \tag{ 4 }$$

where $({\sigma }_{C}^{2},{\sigma }_{V}^{2},{\sigma }_{I}^{2},{\sigma }_{J}^{2},{\sigma }_{H}^{2})$ are the individual photometric uncertainty estimates in the BTS filters and (α, β) the appropriate scale factors for the indices (Equation (1)). Since α² ≫ (1 + β) for these indices (for all populations of interest; Brown et al. 2009), we expect the covariance matrices for most of the stars to generally align with the $[t]$ direction, with only weak uncertainty covariance. Indeed, this is usually the case, though there are exceptions (Figure 18).

Zoom In Zoom Out Reset image size

We are also assuming that the apparent magnitudes and their relevant linear combinations are normally distributed, working in apparent magnitude space rather than flux space because the photometric uncertainties are already reported in magnitudes in the BTS catalog. We impose a photometric uncertainty cut of σ < 0.1 mag (Table 3) to reduce the number of objects that strongly violate this assumption. Nevertheless, long tails in the observed $[t],[m]$ distribution for objects with relatively high photometric uncertainty may be expected.

To estimate the number of components required to best represent the $[t],[m]$ distribution, we employ two commonly used measures, the Akaike Information Criterion (AIC) and the Bayesian Information Criterion (BIC). These measures quantify the badness of fit while penalizing more complex models, with the BIC penalizing overly complex models more severely. More information can be found in Ivezić et al. (2014); these measures take the forms

$$\begin{eqnarray}&&\mathrm{AIC}=2p-2\mathrm{ln}L\end{eqnarray} \tag{ 5 }$$

$$\begin{eqnarray}&&\mathrm{BIC}=p\mathrm{ln}N-2\mathrm{ln}L,\end{eqnarray} \tag{ 6 }$$

where lower values indicate a formally better fit. Here p is the number of parameters in the model, N the number of data points, and L the likelihood (data given model) returned by the mixture modeling procedure. For a GMM consisting of a mixture of K model components representing $q$-dimensional data points, the number of parameters p is given by

$$\begin{eqnarray}&&p=(q\times K)+\left(\displaystyle \frac{q\times K\times (q+1)}{2}\right)+(K-1)\end{eqnarray} \tag{ 7 }$$

so that mixtures with K = 1, 2, 3, 4... model components consist of p = 5, 11, 17, 23... parameters when fitting the 2D $[t],[m]$ distribution. When characterizing the $[m]$ or ($[t],[m]$) distribution with a GMM, we allow K to vary up to large values (usually K = 9) and look for models in which the AIC and BIC stop improving as K is increased.

Figure 19 shows an attempt to reproduce the distribution of $[m]$ only as a GMM (see Appendix B for discussion of the technique). At least two components seem to be required, although the data do not discriminate between the simplest model that fits the data (two components) and a continuum (e.g., eight components). In early trials using data selected only on photometric measurement uncertainty, a mixture model with more than three components would usually include an extremely broad, low-significance Gaussian component. On plotting the $[m]$ counts on a log scale, this component was seen to be fitting handfuls of far outliers in the $[m]$ distribution (with $| [m]| \gt 0.5;$ compare with the range in Figure 19). This may be expected if the outliers are not well represented by the model form; nevertheless, the GMM implementation would attempt to assign a model component to the outliers once the model grew sufficiently complex, which in turn would distort model components much nearer to the location of the main population of objects. Circumventing this outlier problem was the main motivator for outlier removal in $[t],[m]$ when selecting objects for further analysis (Table 3).

Zoom In Zoom Out Reset image size

Figure 20 shows the characterization of the ($[t],[m]$) distribution with a 2D GMM as the number of model components is increased. To examine the impact of changing the number of model components K, the $[t],[m]$ data were split into two equal-size samples (the "training" and "test" sets), and the GMM fit using the "training" set. Samples (of $[t],[m]$) were then drawn from the model and perturbed by measurement covariances ${{\boldsymbol{S}}}_{i}$ from the "test" set, and the ($[t],[m]$) distribution of this predicted set was compared with the "test" set. While models with K = 2, 3, 4 components each provide a reasonable visual match to the observed $[t],[m]$ distribution, the AIC and BIC both indicate that K = 4 provides the best representation of the data, while increasing the number of components beyond K = 4 does not improve the fit further (indeed, the BIC suggests that models with K > 4 fit the data more poorly).

Zoom In Zoom Out Reset image size

Fiber-fed echelle spectroscopy was taken using UVES between 2004 June 22 and 25 (ESO program 073.C-0410(A); PI Dante Minniti). [M/H] estimates were produced in a similar manner to the analysis in Fischer & Valenti (2005) and Valenti & Fischer (2005); typically ∼50 absorption features from a solar spectrum (numerically degraded to the spectral resolution of the observations) are scaled and shifted to find the best match to the observed spectra. In addition to radial velocities, this process also yielded estimates for [M/H] (as well as $\mathrm{log}(g)$ and ${T}_{\mathrm{eff}}$). The [M/H] determination used mainly metal lines, with very few C and O lines in the templates used, which reduces sensitivity in the [M/H] estimates to systematic differences between giants and MS objects (Valenti & Fischer 2005).

The 123 objects in the resulting catalog were trimmed by longitudinal proper motion (μ_l < −2.0 mas yr⁻¹) to produce a sample of 93 likely bulge objects with spectroscopic [M/H] estimates.

Following previous works, which use multicomponent Gaussian mixtures to model the [Fe/H] distributions (e.g., Hill et al. 2011; Schultheis et al. 2017; Zoccali et al. 2017), we also characterize the abundance distribution of the 93 spectroscopically measured likely bulge objects as a Gaussian mixture (Figure 21). Two implementations of GMM with uncertainties are used: the extreme-deconvolution method of Bovy et al. (2011), and scikit-learn XDGMM (Pedregosa et al. 2011). The parameters fitted by the two implementations are generally consistent with each other and are shown in Figure 21 and Table 9.

Zoom In Zoom Out Reset image size

Although the 93 objects have somewhat limited statistical power to distinguish models, it does appear that at least a two-component mixture is preferred. At four or more components, both implementations always include a very broad, almost insignificant component, which suggests overfitting—and indeed the AIC and BIC do not suggest that more than two components are required by these data (Figure 21, right column).

The parameters of the two-component GMM are consistent with those reported by spectroscopic surveys of nearby fields (e.g., Schultheis et al. 2017; Zoccali et al. 2017), both of which find at least two spectroscopic components with similar fractions α_k, centroids, and dispersions. The sample does not include a more metal-poor component that might be suggestive of a halo component (e.g., Ness et al. 2013a; Schultheis et al. 2017).

Additional photometric scatter is simulated as a perturbation in flux. The apparent magnitudes in the "metal-rich" sample are perturbed by amount Δm_p, defined as

$$\begin{eqnarray}{\rm{\Delta }}{m}_{p,i} & = & -2.5{\mathrm{log}}_{10}\left(\displaystyle \frac{{{\rm{F}}}_{0,i}+{\rm{\Delta }}{{\rm{F}}}_{i}}{{{\rm{F}}}_{0,i}}\right)\\ & = & -2.5{\mathrm{log}}_{10}(1+s{ \mathcal N }{(0,1)}_{i}),\end{eqnarray} \tag{ 8 }$$

where ΔF_i is the perturbation in flux, assumed normally distributed, s the scale of the additional flux uncertainty as a multiple of the original unperturbed flux F_0,i, and ${ \mathcal N }{(0,1)}_{i}$ a draw from the unit normal distribution. For large values of s, the normally distributed flux perturbation can cause the perturbed flux values for some simulated objects to go negative; the simulation treats these cases as nondetections and removes affected objects from consideration, thus penalizing simulations with very large simulated flux uncertainty.

Figure 22 shows indications from this test. To aid interpretation in terms of apparent magnitude, we also characterize the sample standard deviation in apparent magnitude caused by the perturbation (which we denote ${s}_{\mathrm{mag}}$), displaying it alongside the input scale s of flux perturbation; the quantity ${s}_{\mathrm{mag}}$ is plotted along the top axes in Figure 22. The rotation curve badness-of-match statistic suggests that observed rotation curve discrepancy can result from increased photometric scatter for scale factor s ≳ 0.35 (in apparent magnitude, ${s}_{\mathrm{mag}}\gtrsim 0.48$), while the $\pi ^{\prime}$ distribution of the "metal-poor" sample is brought into rough agreement with that observed, for scale factor range $0.25\lesssim s\lesssim 0.35$ (or in magnitudes, $0.30\lesssim {s}_{\mathrm{mag}}\lesssim 0.48$).

Zoom In Zoom Out Reset image size

It is difficult to see how the "metal-poor" sample might be subject to such a large additional photometric scatter. For example, the additional photometric scatter is likely far larger than the difference in photometric precision in the two samples from the SWEEPS measurements. Figure 23 shows the internal photometric precision (defined as the rms of the apparent magnitude measurements along the set of images) as a function of apparent magnitude and $\pi ^{\prime}$ for objects in the "metal-rich" and "metal-poor" samples. The "metal-poor" population shows only a slight increase in internal photometric uncertainty compared to the "metal-rich" population, and both are very small (on the order of a few mmag; these objects are well above the photometric completeness limit for the SWEEPS survey). While indeed the internal precision refers to the random component of photometric uncertainty and not the absolute photometric accuracy, a sample difference in photometric uncertainty of ∼0.3–0.5 mag seems highly unlikely for these data.

Zoom In Zoom Out Reset image size

A difference in extinction distribution between the samples, characterized in any way,³¹ if large enough to bring about the ${s}_{\mathrm{mag}}$ ∼ 0.3–0.5 mag additional scatter required, would surely have led to additional observational consequences that are not seen in these data. For example, the observed F814W dispersion of the RCGs in the SWEEPS data set is close to $\sigma ({\rm{F}}814{\rm{W}})\approx 0.17$ mag (Cl08). Even if all this dispersion were due to extinction, which seems unlikely, this would still be a factor ≳2 too low to bring about the observed discrepancies between "metal-rich" and "metal-poor" samples.

When the depth of the bulge along the line of sight is considered, the allowed contribution of differential extinction to $\pi ^{\prime}$ blurring becomes somewhat smaller. For example, assuming that the bulge RCGs are scattered along this line of sight by ±0.50 kpc allows room for only 0.1 mag of photometric blurring due to extinction of any prescription. Since extinction effects would need to apply differentially to the "metal-poor" sample compared to the "metal-rich" sample to bring the two rotation curves into agreement, we conclude that differential extinction effects are likely at least a factor of 3–5 too small to account for the observed rotation curve discrepancy.

It is also not clear why the "metal-poor" sample would be subject to a strongly discrepant extinction distribution (however parameterized) in the first place. The two populations are not strongly different in their projected distributions on the sky (Figure 24), which would seem to argue against, say, the "metal-poor" sample being located within a region on the sky showing stronger, clumpier extinction than the "metal-rich" sample. Additionally, the RCG apparent magnitude distribution in this field does not appear to be bimodal (e.g., Clarkson et al. 2011; Nataf et al. 2013).

Zoom In Zoom Out Reset image size

We point out that this test applies to the dispersion of differential extinction, not to differences in the median extinction between the two samples. Although a difference in median R_V might affect the drawing of the "metal-rich" and "metal-poor" samples using $[m]$, $[t]$ (because those indices are computed in terms of extinction ratios, which are dependent on the prescription for extinction), it would not by itself change the $\pi ^{\prime}$ dispersion for a given population (R_V variations are considered in more detail in Appendix E.5). The $\pi ^{\prime}$ values for "metal-rich" and "metal-poor" populations are both constructed by reference to fiducial ridgelines fit to the observed populations in the SWEEPS CMD. While the interpretation of a given ridgeline with a particular set of population parameters (like [Fe/H], $E(B-V)$, ${{\rm{F}}}_{\mathrm{bin}}$, ${q}_{\min }$, and, to a lesser extent for this population, age) does depend on the median $E(B-V)$, this does not impact the fiducial ridgelines of the observed median populations on the SWEEPS CMD.

We thus reject additional photometric scatter as a cause for the "metal-rich" and the "metal-poor" samples to be drawn from the same kinematic population, because, whatever the cause, its likely magnitude is much too low to have gone unnoticed elsewhere in these data.

If the "metal-poor" sample has a highly discrepant binary fraction or binary companion mass ratio distribution from the "metal-rich" sample, then this might produce a population with larger distance spread, where the additional inferred distance scatter would be biased to closer distances than the mean population—qualitatively similar to the trends observed (e.g., Figure 7).

The binary fraction ${{\rm{F}}}_{\mathrm{bin}}$, minimum binary (initial) mass ratio ${q}_{\min }$, and indeed the shape of the distribution of mass ratio q are not known for the bulge (see, e.g., Calamida et al. 2015) and are difficult to constrain observationally for the sample selected for the present proper-motion study (e.g., Figure 2). A complete search of (${{\rm{F}}}_{\mathrm{bin}}$, ${q}_{\min }$) parameter space, and indeed of the form of the mass ratio distribution, is beyond the scope of the present investigation. Instead, we characterize statistically the distribution of ${\rm{\Delta }}{m}_{\mathrm{bin}}$ due to unresolved binaries, for the CMD region of interest to this study (Figure 2), and draw from this distribution ${\rm{F}}({\rm{\Delta }}{m}_{\mathrm{bin}})$ for each realization of the Monte Carlo trial.

To maximize the impact of a difference in binary fraction between the "metal-rich" and "metal-poor" samples, we assume for the purposes of this test that the "metal-rich" population has no binaries at all, and we perturb it using an unresolved binary fraction to approximate the "metal-poor" population. (This thus allows the excess binary fraction to be tested in the range 0 ≤ ${{\rm{F}}}_{\mathrm{bin}}$ ≤ 1; if we assume that the "metal-rich" sample has a binary fraction of 0.3, then only tests in the range ${{\rm{F}}}_{\mathrm{bin}}$ < 0.7 would be meaningful.) We also assume for this test that it is only the population of unresolved binaries that differs between the two samples (i.e., there is no difference in metallicity distribution between the two samples).

Version 5.0.1 of the BaSTI ³²  suite of simulation tools and stellar population models (Pietrinferni et al. 2004, 2006) is used to produce a representative set of distributions ${\rm{F}}({\rm{\Delta }}{m}_{\mathrm{bin}})$, to characterize for the Monte Carlo draws. Thanks to the capability of BaSTI to accept user-defined random number seeds, simulations that are almost identical except for small changes in input parameters can be run. This allows us to compare synthetic populations on a star-by-star basis, with and without the addition of unresolved binaries.³³

To combine the sophistication of BaSTI with the speed necessary for Monte Carlo trials, the distribution ${\rm{F}}({\rm{\Delta }}{m}_{\mathrm{bin}})$ itself is characterized nonparametrically, using the method outlined in Ivezić et al. (2014, their Section 3.7)—and thus does not depend on a functional form for ${\rm{F}}({\rm{\Delta }}{m}_{\mathrm{bin}})$. This resampling is 10^5–6 times faster than running a BaSTI simulation for each iteration and brings into reach Monte Carlo exploration of the impact of binaries for our purposes here.

BaSTI simulations are run for four choices of the minimum binary initial mass ratio: ${q}_{\min }$ = (0.0, 0.3, 0.5, 0.7). The "bulge" star formation history (Mollá et al. 2000) within BaSTI is used to populate the sample, with scaled-to-solar heavy-element abundances and the Kroupa et al. (1993) initial mass function. Absolute magnitudes are converted to apparent magnitudes using a fiducial distance and reddening. This allows ${\rm{F}}({\rm{\Delta }}{m}_{\mathrm{bin}})$ to be characterized specifically for the population we have selected for proper-motion study. For ${q}_{\min }$ = 0.0, the distribution ${\rm{F}}({\rm{\Delta }}{m}_{\mathrm{bin}})$ turns out to closely resemble ${\rm{F}}({\rm{\Delta }}{m}_{\mathrm{bin}})$ = 1/${\rm{\Delta }}{m}_{\mathrm{bin}}$, while for ${q}_{\min }$ > 0 the distribution becomes more complicated and nonparametric resampling is preferred (Figure 25).

Zoom In Zoom Out Reset image size

In none of the cases (${q}_{\min }$ = 0.0, 0.3, 0.5, 0.7) do we find that the presence of an additional binary population can account for the difference between the observed "metal-rich" and "metal-poor" rotation curves. (Figure 26 shows the cases ${q}_{\min }$ = 0.0 and ${q}_{\min }$ = 0.7). Only the skewness of the $\pi ^{\prime}$ distribution ever approximates that of the "metal-poor" population (at ${{\rm{F}}}_{\mathrm{bin}}$ ≳ 0.5), while the rotation curve and $\pi ^{\prime}$ spread do not overlap for any binary fraction.

Zoom In Zoom Out Reset image size

We therefore conclude that an excess of unresolved binaries in the "metal-poor" over the "metal-rich" population is highly unlikely to be responsible for the difference in rotation curves.

To estimate the differential scatter in photometric parallax produced by differing [Fe/H] dispersions between "metal-rich" and "metal-poor" samples, a synthetic composite stellar population is produced for the SWEEPS field by sampling BaSTI simulations (computed for all three cameras and resampled in the manner of Appendix D.2), which include the effects of age, [Fe/H] spread, and unresolved stellar binaries. The synthetic populations are perturbed by photometric uncertainty (in all seven filters), photometric parallax, and reddening, where the width of the distributions in all three quantities can be specified separately for each population.

This produces a SWEEPS CMD and $[t],[m]$ distribution for the synthetic population. Synthetic objects are selected for further "study" in a similar manner to that for the real data (e.g., Table 3); in particular, synthetic SWEEPS CMD objects must fall within the selection box in the SWEEPS filters (Figure 2). The surviving synthetic objects are then classified as likely "metal-poor" and "metal-rich" populations in the same manner as for the observed data (using the GMM components in $[t],[m]$ that were fitted to the real data), isolating "observed" samples of "metal-rich" and "metal-poor" objects. In this manner, the synthetic samples are isolated in a similar fashion to those drawn from the real data.

Finally, best-fit loci are determined for the model absolute magnitudes of the synthetic "metal-rich" and "metal-poor" samples, and the differences ΔM_V from these loci are determined for every object in the samples. The model absolute magnitude is used rather than the apparent magnitude because we wish to isolate the impact of metallicity spread on intrinsic magnitude scatter—i.e., before distance, reddening, and photometric uncertainty have perturbed the measurements (which impacts the sample selection), but including the intrinsic effects of age, [Fe/H], and binarity.

Modeling the selection cuts on the synthetic samples requires simulating the composite stellar population of the SWEEPS field. This field is somewhat complex, consisting of at least three distinct populations (bulge, local disk, halo), each of which could well consist of multiple subpopulations or a continuum.

Full population decomposition presents a formidable challenge (e.g., Gennaro et al. 2015) and is complicated by the difficulty in adequately accounting for extinction across the broad wavelength range of the BTS photometry in the inner bulge region (e.g., Nataf et al. 2016). To produce a reasonable approximation to the selection effects at work in the SWEEPS field, a multicomponent stellar population is instead simulated with parameters drawn from the literature and the [Fe/H] spread estimated in this work (Appendix C). A set of about a dozen synthetic populations with various parameter settings are simulated using BaSTI, with typically five components from this set combined appropriately to produce a synthetic composite population for the SWEEPS field, with mixture parameters tuned by hand to provide an approximate match to the observed SWEEPS CMD and $[t],[m]$ distribution.

All population components used the same prescription for binaries, with binary fraction 0.35 and minimum binary mass ratio 0.0. The initial mass function followed the Kroupa et al. (1993) prescription for all components over the BaSTI default mass range ($0.1\leqslant M/{M}_{\odot }\leqslant 120$). Convective core overshooting was not selected for any model component, and mass-loss parameter η = 0.4 was used throughout. When not using a predetermined star formation history supplied by BaSTI, the star formation histories were specified as a series of single bursts at given ages, with [Fe/H] described as a Gaussian with a user-specified centroid and standard deviation. Specific details for various population components follow below.

For the foreground disk, the formation history of Rocha-Pinto et al. (2000) was used (the default "local disk" scenario within BaSTI), typically forming 5%–10% of the stars in the simulation sets.

Stellar halo components were simulated using the bimodal [Fe/H] distribution reported by An et al. (2013) from Sloan Digital Sky Survey photometry; this model consists of a very metal-poor component centered at [Fe/H] ≈ −2.33 and another slightly less metal-poor component centered at [Fe/H] ≈ −1.67. For a bimodal bulge population following any of the GMM fits to our spectroscopic data, or for the [Fe/H] distribution of Zoccali et al. (2017) near the SWEEPS field, this separate halo component is necessary to populate the regions in $[t],[m]$ space for objects with [Fe/H] ≲ −2.0.

Bulge components were constructed separately as normally distributed [Fe/H] distributions specified through the BaSTI web interface, using the characterization presented in Appendix C, for both the two- and three-component GMM decompositions. For components less metal-rich than [Fe/H]₀ < + 0.3, separate runs were simulated using the "scaled-to-solar" and "α-enhanced" options within BaSTI in order to allow some exploration of α-enhancement on population spread in the $[t],[m]$ diagram. A variety of age prescriptions were attempted, mostly to improve the match at the bright end of the SWEEPS CMD, either by ascribing a single burst of star formation to each metallicity or by assigning several bursts to each metallicity (e.g., bursts at 5.0, 6.0, and 7.0 Gyr for a component with [Fe/H]₀ = −0.42). We have not yet explored more sophisticated age–metallicity prescriptions through user-defined star formation histories (e.g., Haywood et al. 2016; Bensby et al. 2017).

The more continuous bulge star formation history of Mollá et al. (2000, used as a default in BaSTI) was also tried, for "scaled-to-solar" isochrones, for "α-enhanced" isochrones, and for varying admixtures of the two.

We have not yet explored a separate "thick-disk" component in this context. The metal-poor wing of the bulge distribution or the metal-rich wing of the halo component could mimic such a population in the $[t],[m]$ diagram, and we do not make the distinction here.

Figures 27 and 28 show examples of the synthetic populations thus produced. None of the population mixtures that we have produced quite reproduce both the observed SWEEPS CMD and the $[t],[m]$ diagram, although in view of both the challenges of extinction characterization and apparent simulation truncations imposed by BaSTI itself (Appendix E.3), full reproduction is likely to be difficult. The basic two-component bulge we simulate here produces an $[t],[m]$ distribution that is much more strongly bimodal than that observed (e.g., Figure 3), while the 10-component "bulge" star formation history within BaSTI (Mollá et al. 2000) produces an $[t],[m]$ distribution that is too smooth compared to that observed.

Zoom In Zoom Out Reset image size

Several methods were attempted to bring the simulated $[t],[m]$ distribution into closer agreement with that of the observed data in Figure 3. One simple ansatz is to simply multiply the BTS estimated photometric uncertainties by a factor of 2 before selection and computation of $[t],[m]$ (Figure 28, top right panel). Another is to apply Gaussian blurring in $[t]$ and $[m]$ separately (bottom middle panel of Figure 28). Varying R_V with a Gaussian of width ${\sigma }_{{R}_{V}}=0.52$ does bring the marginal distribution reasonably close to that observed (bottom right panel of Figure 28), although the $[t],[m]$ distribution that results is distorted compared to the observed sample (particularly the "metal-rich" sample), and in addition the required ${\sigma }_{{R}_{V}}$ is at least a factor of ∼2 larger than that suggested by the SWEEPS CMD (Appendix E.5, which also shows the $[t],[m]$-blurring effect due to R_V variations that are compatible with the SWEEPS data).

Zoom In Zoom Out Reset image size

For the purposes of estimating the impact of varying [Fe/H] distribution on relative photometric parallax variations, we retain the two-component bulge model with and without BTS uncertainty scaling, for further investigation; the former is consistent with estimated [Fe/H] distributions and estimates of photometric uncertainty, while the latter is the "broadened" option among those tried that closely resembles the observed distribution (Figure 3).

While conducting tests on the simulated data sets, it quickly became apparent that samples generated with the current version of BaSTI³⁴ show truncation at extremes of both high and low metallicity, leading to a hard edge in the CMD of the simulated population that has no counterpart in the reported [Fe/H] distribution. This truncation, characterized in Appendix F, impacts the metal-rich simulated bulge sample more strongly than its metal-poor simulated counterpart and thus could artificially enhance the discrepancy in absolute magnitude breadth between the metal-rich and metal-poor simulated components.

This hidden systematic complicates efforts to characterize the excess magnitude scatter due to differing [Fe/H] distributions, with much of the most metal-rich end of the metal-rich simulated sample assigned apparently incorrect magnitudes (absolute and apparent). We therefore adopt a restricted-sample estimate of the magnitude scatter, by sampling only the fainter side of the magnitude distribution for both samples in the comparison. Specifically, we use the quantity σ_hi defined by³⁵

$$\begin{eqnarray}&&{\sigma }_{\mathrm{hi}}^{2}\equiv \displaystyle \frac{1}{N(m\geqslant \overline{m})}\displaystyle \sum _{m\geqslant \overline{m}}{({m}_{i}-\overline{m})}^{2},\end{eqnarray} \tag{ 9 }$$

where, for the special case of a large, strictly symmetric distribution, σ_hi closely approximates the sample standard deviation. A practical challenge is to identify the median magnitude $\overline{m}$ from a truncated asymmetric distribution. For these simulations, $\overline{m}$ is estimated by discarding the most negative Δm samples (thus discarding objects near and outside the truncation limits) and fitting a Gaussian function to the histogram of Δm values. This fit is only used to estimate $\overline{m}$, which thus allows σ_hi to be estimated following Equation (9). This then allows the restricted-sample scatter σ_hi to be estimated for the metal-poor and metal-rich samples separately, and the excess difference to be characterized as the quadrature difference between the two.

The final step is then to convert the excess scatter ${\sigma }_{\mathrm{hi}}$ estimated from the simulated population components to the additional flux scatter s felt by the metal-poor sample compared to the metal-rich sample. To enable this conversion, the relationship between restricted-sample scatter σ_hi and the flux perturbation scale s that generated it was determined by simulation. Synthetic populations with perturbation flux distribution were produced following Equation (8), subject to the same censoring for negative flux as before (Appendix D.1). The apparent magnitude scatter σ_hi was then found for each synthetic population as described above (and as performed for the simulated BaSTI data sets). Finally, the relationship between ${\sigma }_{\mathrm{hi}}$ and s was characterized by fitting a seventh-order polynomial in both directions. Table 10 and Figure 29 show this characterization. This allows us to relate the restricted-sample scatter found from BaSTI simulations back to the flux ratio perturbation scale s, and finally to compare the scale of the perturbation suggested by differing [Fe/H] distributions to the additional scale of flux perturbations s that our observational data would require if the "metal-poor" sample really were a blurred version of the "metal-rich" sample.

Zoom In Zoom Out Reset image size

We are finally in a position to estimate the additional scatter in absolute magnitude due to differential metallicity scatter. Figure 30 shows the results of applying the selection criteria to the simulation including the two-component bulge model, a two-component halo, and the local disk component.

Zoom In Zoom Out Reset image size

Figure 31 illustrates the characterization of absolute magnitude scatter σ_hi, while Table 11 shows the evaluation of the excess flux scatter s for "metal-poor" compared to "metal-rich" samples. Two simulated populations were evaluated in this manner: one including the two-component bulge model, and the other with the BTS uncertainties multiplied by a factor of 2 before selection to broaden the distribution in $[t],[m]$. In both cases, the excess fractional flux scatter s is less than 0.1; we find s ≈ 0.09 for the two-component bulge model, while s ≈ 0.07 for the enhanced-uncertainty version of this model.

Zoom In Zoom Out Reset image size

We contacted the authors of the BaSTI web tools regarding its internal truncation (detailed in Appendix F). In response, S. Cassisi (2017, private communication) kindly added a high-metallicity point to BaSTI's internal metallicity grid (since in BaSTI version 5.0.1 the metallicity range covered by the simulator is more restrictive than that covered by the isochrone set) and recomputed sets of synthetic populations using the updated version of the simulator.³⁶ Visual inspection of the $[t],[m]$ distribution and the SWEEPS CMD drawn from the Cassisi simulations indicates similar behavior to those from v5.0.1, except without the sharp edges truncating the metal-rich end of the synthetic population.

In this paper we retain the statistics derived using BaSTI v5.0.1 since that is the version currently available to the community. However, the comparison with the Cassisi version is instructive. Application of the half-sample techniques of Appendix E.3 to both the Cassisi and v5.0.1 simulations yielded highly similar results (${\sigma }_{\mathrm{hi}}$ differing by <4%), as might be expected since this measure uses the side of the ΔM distibution far from the truncation limit. The Cassisi simulations also allow a direct estimate of the accuracy of the one-sided measure adopted in Appendix E.3, by comparing ${\sigma }_{\mathrm{hi}}$ to the ΔM standard deviation of the objects in the dominant component of the Cassisi simulation (see Figure 30 for the dominant and "background" components for metal-rich and metal-poor simulated populations). In the Cassisi simulations, the ΔM standard deviation is roughly 20% smaller than the estimate σ_hi, suggesting that our estimates of the excess photometric scatter in Table 11 may be overestimates.

We therefore find that the combination of differing metallicity spreads between "metal-poor" and "metal-rich" samples and differing selection effects in both the $[t],[m]$ distribution and SWEEPS CMD contributes differential flux scatter that is not larger than σ_hi ≈ 0.1 mag, or additional flux standard deviation s ≈ 0.08. This additional scatter is a factor of 3 too small to bring the observed "metal-poor" and "metal-rich" proper-motion-based rotation curves into agreement by itself (Figure 22), and we conclude that the apparent difference in proper-motion rotation curves between the two samples is not an artifact of differences in the underlying [Fe/H] distribution.

As a second check, we can compare the [Fe/H] distribution of the objects classified as "metal-poor" and "metal-rich" with the simulated [Fe/H] values for the relevant bulge model components. We find that indeed the misclassification rate in this synthetic population-based simulation appears to be low (Figure 32). Possible contaimination is explored further in a purely empirical manner in Appendix G.

Zoom In Zoom Out Reset image size

The framework of this appendix also allows us to investigate the impact of R_V variations on $[t],[m]$-based determinations. The extinction-free indices $[t],[m]$ assume a particular extinction prescription (Cardelli et al. 1989, using ${R}_{V}=2.5$). While $[t],[m]$ are therefore insensitive to variations in $E(B-V)$ for a particular value of R_V, variations in ${R}_{V}$ could impact the distribution of points in the $[t],[m]$ diagram, by altering the relationships between apparent magnitudes in the BTS filters and those assumed when computing $[t],[m]$.

We appeal to the SWEEPS CMD to estimate limits on the magnitude of R_V variations in this field. Assuming that the distance distribution due to the physical depth of the bulge can in this field be characterized by a Gaussian with width parameter σ_d kpc, the observed apparent magnitude scatter of RCGs in this field then sets an upper limit on ${R}_{V}$ variations for assumed $E(B-V)$. In the SWEEPS data set, the observed F814W dispersion of the RCG is $\sigma ({\rm{F}}814{\rm{W}})\approx 0.17$ mag (Cl08).

For this appendix we adopt $E(B-V)$ = 0.5 (Ca14) as a representative value (the implied R_V variations would become smaller for larger $E(B-V)$). The extreme case of distance dispersion, σ_d = 0, then admits ${R}_{V}$ variation of ${\sigma }_{{R}_{V}}\approx 0.45$. However, the bulge has nonzero depth along the line of sight; picking a representative distance distribution of ${\sigma }_{d}\approx 0.5\,\mathrm{kpc}$ suggests that variation closer to ${\sigma }_{{R}_{V}}$ ≈ 0.25 is more likely. Both estimates for ${\sigma }_{{R}_{V}}$ are conservative upper limits, since they ascribe none of the observed RCG apparent magnitude dispersion to photometric uncertainty, luminosity variations within the RCG sample, or $E(B-V)$ variation.

To estimate the impact of ${R}_{V}$ variation on the $[t],[m]$ distribution (and thus sample selection and cross-contamination), a synthetic population was constructed using BaSTI population components tuned to the estimated metallicity distribution for this field. Full details of this procedure, which was implemented to explore metallicity-dependent selection and characterization systematics (Section 5.1), can be found in Appendix E.

Figure 33 shows the comparison of a simulated $[t],[m]$ population, with and without ${R}_{V}$ variations at the ${\sigma }_{{R}_{V}}=0.25$ level admitted by the SWEEPS data set. For each relevant WFC3 filter, the scale factors ${A}_{X}/E(B-V)$ were estimated by linear interpolation in ${R}_{V}$ using information shown in Table 4. The simulated magnitudes were thus perturbed into "observed" magnitudes using different ${R}_{V}$ values for each star, but the $[t],[m]$ were computed using the $\alpha ,\beta$ values appropriate for ${R}_{V}=2.5$. This then mimics the use of a single ${R}_{V}$ value to compute $[t],[m]$ for a population that in reality shows ${R}_{V}$ variations.

Zoom In Zoom Out Reset image size

Comparing the synthetic $[t],[m]$ distributions with and without ${R}_{V}$ variations (Figure 33), it seems unlikely that ${R}_{V}$ variations at the level admitted by the SWEEPS CMD can contribute a strong effect on GMM fitting or sample selection in $[t],[m]$; the impact of ${R}_{V}$ variations is simply too small. We therefore proceed under the assumption that indeed ${R}_{V}\approx 2.5$ for all objects in the SWEEPS field of view.

To investigate whether BaSTI is applying the truncation to the median population in a requested sample, test populations were simulated for bursts of star formation of equal magnitude but with very narrow [Fe/H] distributions. Figure 34 shows an example for a scaled-to-solar set of isochrones, with [Fe/H] = {−3.0, −2.5, −2.0, −1.5, −1.0, −0.5, +0.0, +0, 4, +0.5}, all with spread σ_[Fe/H] = 0.0001 dex to isolate selection effects applied to the mean populations in each case. The two most metal-poor populations and the single most metal-rich population are found to be forced away from their specified values, probably to some internal limit. Reading off the figure, the most metal-poor populations seem to be brought up to [Fe/H] ≈ −2.3, with the most metal-rich brought down to [Fe/H] ≈ +0.40. These values are entirely consistent with the [Fe/H] limits suggested by the BaSTI documentation referenced earlier.

Zoom In Zoom Out Reset image size

This suggests that BaSTI enforces [Fe/H] limits on the median populations requested in a simulation.

To investigate whether BaSTI applies a truncation to [Fe/H] values that are carried outside internal [Fe/H] limits owing to the specified population spread, test populations were simulated including a single population well away from the limits and one component each just inside the two limits. Components were specified with $[\mathrm{Fe}/{\rm{H}}]=\{-2.5,-1.2,+0.05\}$, all with specified spread σ_[Fe/H] = 0.1 dex, to ensure that the two components near the [Fe/H] limits each include substantial numbers of objects outside these limits, while the middle population has very few such objects.

Figure 35 shows the resulting simulation. Curiously, although the [Fe/H] values reported in the simulated populations show no truncation, the CMD and the simulated absolute magnitudes quite clearly do show truncation, with a hard edge at both the upper and lower [Fe/H] extrema.

Zoom In Zoom Out Reset image size

We therefore find that BaSTI does not truncate [Fe/H] values at the stage of assignment to simulated objects, and these nontruncated [Fe/H] values are carried through to the output simulated population. However, a truncation is applied at some stage before the absolute magnitudes are included in the simulated population. This results in both a hard edge to the distribution of simulated absolute magnitudes and a discrepancy between the reported [Fe/H] values and the absolute magnitudes in the simulation output.

To chart the behavior of the truncation near the [Fe/H] limits in more detail, we simulated a single test population near the metal-rich limit. Figure 36 shows the result for a scaled-to-solar component with [Fe/H] = +0.24 and scatter σ_[Fe/H] = 0.19. In this case, the truncation appears to be quite dramatic, with a narrow, highly overrepresented component in the ΔM distribution.

Zoom In Zoom Out Reset image size

However, the behavior of the simulator near an [Fe/H] limit is not as straightforward as a simple substitution of the [Fe/H] limit for all objects beyond it. Figure 31 shows a simulated metal-rich population partitioned by [Fe/H], which allows us to distinguish objects that were assigned [Fe/H] values above the metal-rich limit (and thus would be assumed to be truncated). Objects with outlier [Fe/H] values do not only appear at the location where absolute magnitudes pile up; a substantial fraction show magnitudes deeper into the main population (see the bottom left panel of Figure 31).

That the pileup implying truncation is also observed at the metal-rich edge of the population with simulated metallicities within the limits according to the BaSTI documentation suggests that the effective metallicity limits may differ from those documented; see the middle left panel of Figure 31. (We have not yet dissected this simulated population by binarity, which might offer another avenue for objects to wander into truncation territory.)

We therefore find that the internal truncation applied by BaSTI is not limited to a simple pegging of values to an internal boundary. The behavior probably necessitates some sort of selection on ΔM to produce a cleaner unaffected sample. We adopt one such approach in Appendix E.3.

The BTSv1 photometric catalog was produced using daophotII on summed images in each filter before combining into the final catalog (Section 2.2 and references therein), while the BTSv2 catalog uses "effective PSF" methods (e.g., Anderson & King 2006), the details of which vary depending on the brightness regime of the object used. Objects in the brightness range of our proper-motion sample (Figure 2) were measured using the kstwo code by J. Anderson, which fits position and flux for each star across all exposures simultaneously (see Bellini et al. 2017, for details). For sufficiently bright and isolated objects, the source position and flux were fit independently, while for fainter and/or less isolated stars the flux was measured using forced photometry (with the star position fixed). Because both the source brightness and degree of isolation depend on the filter and camera used, a given star might be measured using forced photometry in some filters but not others.

Figure 38 presents the comparison of apparent magnitude between BTSv2 and BTSv1, over the apparent magnitude range of interest to the present work (stars were cross-matched between the two catalogs using their equatorial coordinates). The random component of the apparent magnitude difference is $\lesssim 0.06$ mag in all filters for most objects (compare with the apparent magnitude selection criteria in Table 3). Systematic offsets between the data sets are less than 0.02 mag.

Zoom In Zoom Out Reset image size

Figure 39 presents the star-by-star proper-motion comparison between the BTSv2 and SWEEPS catalog for objects in our apparent magnitude range of interest. Any difference in scale between the proper-motion determinations is below 1%. A small offset ${\rm{\Delta }}{{\boldsymbol{\mu }}}_{{\rm{F}},0}$ ≈ $(0.3,0.1)$ between the two catalogs is apparent, as expected if the proper-motion zero-point of the two catalogs depends ultimately on the differing depth of the two surveys. The proper-motion differences show rms scatter $\approx 0.3$ mas yr⁻¹. Assuming that the full SWEEPS proper motions carry uncertainty ${\epsilon }_{\mathrm{SWEEPS}}\lesssim 0.12$ mas yr⁻¹ (Appendix A), the BTSv2 proper motions for this field would contribute approximately _BTS ≈ 0.27 mas yr⁻¹. No trend in the proper-motion differences was found against apparent magnitude, proper motion, or position.

Zoom In Zoom Out Reset image size

To check the orientation of the astrometric frames of the catalogs, positions in the SWEEPS and BTSv2 catalogs were matched to their entries in the First Gaia Data Release, which should provide absolute positions on the International Celestial Reference System (ICRS) in the 2015.0 epoch, albeit possibly with residual distortions in these crowded regions (Gaia Collaboration et al. 2016a, 2016b). Objects in the Gaia apparent magnitude range $18.0\leqslant G\leqslant 19.5$ were selected for cross-matching, as a trade-off between quality of Gaia measurement and the desire to avoid highly saturated objects in the SWEEPS and BTS catalogs; this leaves a few thousand objects with which to probe positional differences. To minimize random scatter in the comparison, positions from the SWEEPS and BTSv2 catalogs were advanced to their positions in the 2015.0 epoch using the measured proper motions in each catalog; the Gaia DR1 catalog does not contain proper motions for these objects.

Figure 40 maps the astrometric offsets from the Gaia DR1 frame for both the SWEEPS and BTS catalogs. While the two catalogs are slightly offset with respect to the Gaia DR1 frame (by $\lesssim 0\buildrel{\prime\prime}\over{.} 3$ in each coordinate), no rotational flow pattern is detected that would suggest misalignment of either of the reference frames. The scale of the positional residuals is surprisingly large, with flow pattern common to both catalog comparisons that reaches up to $\sim 0\buildrel{\prime\prime}\over{.} 15$ in some regions. The largest residual structure is found at a similar location in the comparisons to both the SWEEPS and BTSv2 catalogs, despite the two HST catalogs being taken at different camera orientations and with different field centers, so we suspect that the flow pattern is dominated by distortion in the Gaia DR1 frame in these crowded regions. (Comparison of Subaru measurements with Gaia DR1 positions near the core of the Sextans dwarf galaxy also shows a flow pattern of offsets on a scale of $\sim 50^{\prime\prime}$ with a gap in Gaia DR1 coverage; C. Dinescu 2018, private communication.) We expect that this flow pattern will vanish in comparisons to future Gaia data releases that have included the more sophisticated treatment for crowding outlined in Pancino et al. (2017). Based on the scale of the flow patterns near the corners of the difference maps (≲01), we conclude that the astrometric reference frames of the SWEEPS and BTSv2 catalogs are aligned with the ICRS frame to better than ∼005.

Zoom In Zoom Out Reset image size

Finally, to investigate whether our results qualitatively change when moving from BTSv1 to BTSv2, we have performed a preliminary re-analysis using the BTSv2 measurements only, following the procedures of Section 3 as far as the production of the proper-motion rotation curves.

Not all the selection steps are common to both catalogs; for example, BTSv2 does not contain F606W measurements (as were used when we melded the SWEEPS and BTSv2 catalogs in Section 3), which thus alters the initial selection of objects, and BTSv2 contains additional information that can be used to select objects by measurement quality (details can be found in the BTSv2 README file). Additionally, BTSv2 proper-motion uncertainties have not yet been characterized as fully as the ACS/WFC uncertainties in the SWEEPS field (e.g., Calamida et al. 2015). Finally, the proper-motion zero-points of the two catalogs differ, and the appropriate value of ${D}_{0}$ to use for BTSv2 has not yet been established.

Figure 41 shows the results. While some of the fine structure in the rotation curves appears to differ when compared to the BTSv1-based analysis, the behavior we observe is not substantially changed by use of the BTSv2 catalog; the "metal-rich" and "metal-poor" rotation curves still differ, with the "metal-rich" curve showing a steeper gradient. Full development of the chemically dissected bulge rotation curves will be reported in a future communication after the work has been extended to all four BTS fields.

Zoom In Zoom Out Reset image size