ARGOS – III. Stellar populations in the Galactic bulge of the Milky Way

Abstract

We present the metallicity results from the ARGOS spectroscopic survey of the Galactic bulge. Our aim is to understand the formation of the Galactic bulge: did it form via mergers, as expected from Λ cold dark matter theory, or from disc instabilities, as suggested by its boxy/peanut shape, or both? Our stars are mostly red clump giants, which have a well-defined absolute magnitude from which distances can be determined. We have obtained spectra for 28 000 stars at a spectral resolution of R = 11 000. From these spectra, we have determined stellar parameters and distances to an accuracy of <1.5 kpc. The stars in the inner Galaxy span a large range in [Fe/H], −2.8 ≤ [Fe/H] ≤ +0.6. From the spatial distribution of the red clump stars as a function of [Fe/H], we propose that the stars with [Fe/H] > −0.5 are part of the boxy/peanut bar/bulge. We associate the lower metallicity stars ([Fe/H] < −0.5) with the thick disc, which may be puffed up in the inner region, and with the inner regions of the metal-weak thick disc and inner halo. For the bulge stars with [Fe/H] > −0.5, we find two discrete populations: (i) stars with [Fe/H] ≈ −0.25 which provide a roughly constant fraction of the stars in the latitude interval b = −5° to −10°, and (ii) a kinematically colder, more metal-rich population with mean [Fe/H] ≈ +0.15 which is more prominent closer to the plane. The changing ratio of these components with latitude appears as a vertical abundance gradient of the bulge. We attribute both of these bulge components to instability-driven bar/bulge formation from the thin disc. We associate the thicker component with the stars of the early less metal-rich thin disc, and associate the more metal-rich population concentrated to the plane with the colder more metal-rich stars of the early thin disc, similar to the colder and younger more metal-rich stars seen in the thin disc in the solar neighbourhood today. We do not exclude a weak underlying classical merger-generated bulge component, but see no obvious kinematic association of any of our bulge stars with such a classical bulge component. The clear spatial and kinematic separation of the two bulge populations (i) and (ii) makes it unlikely that any significant merger event could have affected the inner regions of the Galaxy since the time when the bulge-forming instabilities occurred.

1 INTRODUCTION

Recent observations show that the stars in the Galactic bulge region comprise several kinematic and chemical components spanning a range in α-enhancement (Babusiaux et al. 2010; Bensby et al. 2010). Although the bulge region can no longer be described as a single stellar population, it appears to be a mostly old, α-enhanced system whose stars formed over a short period of time early in the life of the Galaxy (Zoccali et al. 2003; Lecureur et al. 2007; Meléndez et al. 2008; Alves-Brito et al. 2010; McWilliam & Zoccali 2010).

Near-infrared (NIR) imaging (Okuda et al. 1977) and COBE/Diffuse Infrared Background Experiment (Dwek et al. 1995; Smith, Price & Baker 2004) revealed the boxy/peanut morphology of the Milky Way bulge. Blitz & Spergel (1991) established the bar-like structure of the inner bulge region that is now understood to go with the boxy/peanut bulge structure. The major axis of this bar points into the first Galactic quadrant, at about 20° to the Sun–centre line. Its semilength is about 3.1–3.5 kpc (Gerhard 2002) and its axial ratios are about 1: 0.33: 0.23 (Dwek et al. 1995). NIR star counts (e.g. López-Corredoira, Cabrera-Lavers & Gerhard 2005) provide evidence for a longer flatter component of the bar, but the nature of this long bar is uncertain. Martinez-Valpuesta & Gerhard (2011) and Romero-Gómez et al. (2011) argue that the apparent long bar is an artefact associated with leading spiral features at the end of the shorter primary bar.

In Λ cold dark matter (ΛCDM) simulations of Galaxy formation (e.g. Abadi et al. 2003; Kobayashi & Nakasato 2011), bulges are built up primarily through hierarchical mergers. However, it is now widely believed from simulations and observations that the boxy/peanut-shaped bulges, like the bulge of the Milky Way, are not merger products but formed via instability of the inner disc (see e.g. Combes & Sanders 1981; Raha et al. 1991; Bureau & Freeman 1999). The early simulations showed that after a few revolutions, the flat disc develops a bar which then becomes vertically unstable and puffs up into a boxy/peanut bulge structure. More recent N-body simulations of bulge formation via disc instability (Athanassoula 2008) help to understand the time-scales and the spatial, kinematic and chemical abundance tracers of this instability process. For example, the simulations can be used to interpret the mapping of the early disc into the boxy/peanut bulge via the instabilities.

The simulations of bar/bulges produced by disc instability also show that cylindrical rotation is expected, and these are supported by observations of edge-on boxy/peanut bulges (e.g. Falcón-Barroso et al. 2006). If the Galactic boxy/peanut bulge is a disc instability product, then we would expect it to show cylindrical rotation. From their survey of M giants in the inner Galaxy and based on such kinematical criteria, Howard et al. (2009) constrain any merger-generated component of the bulge to be less than 8 per cent. Several authors continue to argue in favour of the importance of dissipational collapse or mergers in the formation of the Milky Way bulge (Zoccali et al. 2008b; Babusiaux et al. 2010), based on their interpretation of new chemical and kinematic data.

We have completed a large spectroscopic survey of about 28 000 stars in the bulge and surrounding disc in order to measure kinematics and chemistry and so to constrain the formation process of the bulge. We aim to establish whether it is primarily a merger or an instability product. We used models of classical bulges to estimate that 28 000 stars are sufficient to detect a 5 per cent merger-generated bulge underlying an instability-generated bar/bulge. We note that this estimate is contingent upon the classical bulge being a slowly rotating and dispersion-supported component. We would not be able to identify a classical bulge component which has been spun up into cylindrical rotation as in the simulations of Saha, Martinez-Valpuesta & Gerhard (2012).

This is the third paper in the ARGOS series and concentrates primarily on the metallicity distribution function (MDF) for the ARGOS stars in the inner Galaxy. Paper I (Ness et al. 2012) examined the nature of the bimodal distribution of apparent magnitudes of the red clump stars along the minor axis of the bulge. This split clump appears to be a generic feature of boxy/peanut bulges and provides a useful guide to the interpretation of the metallicity components in the bulge. Paper II (Freeman et al. 2013) gave a description of the survey, our choice of fields, and the details of our observations, sample selection and analysis, including estimation of stellar parameters and distances. Paper IV (Ness et al., in preparation) will describe the kinematical properties of the bulge as evaluated from the ARGOS stars. In Paper V, we will further interpret the origin of the MDF components by using our dynamical models to discuss how the instability of the disc generates a phase-space-dependent mapping of the disc into the bulge.

In Section 2, we examine the spatial and metallicity distributions of stars in our sample, identify the components in the metallicity distribution in the inner Galaxy, and discuss the abundance gradients of these components in latitude and galactocentric cylindrical radius R_G. In Section 3, we present the components that comprise our MDF and discuss their properties in Section 4. We present our results for the α-enhancement of the bulge as a function of latitude in Section 5, and show the relationship between α-enhancement and [Fe/H] as a function of R_G. In Sections 6, 7 and 8 we discuss the interpretation and origin of the MDF components that we have identified and compare our results to other studies. Section 9 gives our conclusions.

2 RESULTS

2.1 Overview

For reference, Fig. 1 shows the location of the 28 two-degree fields of the ARGOS survey. Each field has about 1000 stars, mostly red clump giants, as discussed in detail in Freeman et al. (2013).

As part of the survey, we have obtained radial velocities for each star with an error of <0.9 km s⁻¹ and determined the stellar parameters T_eff, log g, [Fe/H] and [α/Fe] for each of our stars. The [Fe/H] values are determined from the Fe lines in the Ca-triplet region of the spectrum, using a masked χ² technique against a grid of model spectra. We have excluded stars with poor quality or low-signal-to-noise-ratio (S/N) spectra, giving a final sample of about 25 500 stars. The resolution of the spectra is about 11 000 and the typical S/N for our stars is 75 per resolution element. From calibration against stars in well-observed globular and open clusters, the typical accuracies of the [Fe/H] and log g values are about 0.09 and 0.30 dex, respectively. Fig. 2 show examples of the (T_eff, log g) distribution of stars for our minor axis field at (l, b) = (0°, −10°) with BaSTI isochrones (Cassisi et al. 2006) of age 10 Gyr superimposed. The six panels show the distributions and isochrones for six intervals of [Fe/H] from −2.5 to +1.0. The isochrones have the mean abundance for the panel range in [Fe/H], and [α/Fe] = 0.4 for the isochrones at [Fe/H] < −0.25. The scatter of the stars about the isochrones is as expected. For some stars, the stellar parameters returned by our χ² process were too close to the edge of the model grid (i.e. log g < 0.5) and these stars were excluded from the analysis.

In order not to exclude metal-poor giants from the survey, we selected stars with (J − K)_○ > 0.38 and therefore expect to find some nearby dwarfs with spectral types later than about G2, as seen in Fig. 2. These dwarfs are readily identified from their position relative to the isochrones in the (T_eff, log g) plane. We see that most of the stars are clustered around the red clump (log g = 2–3 and T_eff = 4500–5300 K), as we would expect (Zhao, Qiu & Mao 2001; Valentini & Munari 2010; Saguner et al. 2011). Zhao et al. (2001) studied Hipparcos clump giants at high resolution and found two populations: a metal-rich clump for which the gravity is higher ([Fe/H] > 0, log g = 2.5–3.1) plus a few more metal-poor stars with lower gravity ([Fe/H] < −0.3, log g = 1.9–2.4). We also find that the gravity of the main clump of stars is shifted towards lower values for lower [Fe/H], and we expect that the clump is not contributing significantly to the stellar population for [Fe/H] < −1.0.

The clustering around the red clump region is seen across all of the fields in our survey. Fig. 3 shows the density distribution in the (T_eff, log g) plane for the 25 500 stars in our survey, plotted in bins of metallicity against the isochrones. While the more metal-poor stars lie close to the isochrones, the more metal-rich stars deviate progressively more and more from the isochrones. The deviations are larger for the giants that lie above the red clump and are most marked for stars in the most metal-rich bin, [Fe/H] > 0.25. We need to discuss these observed deviations from the isochrones, because we will be deriving isochrone distances for all of the survey stars. The apparent deviations of our stars to the left may be due to several potential sources of error. Before discussing these error sources, we briefly summarize the relevant details of our procedure for deriving stellar parameters (see Freeman et al. 2013 for further discussion).

While our χ² technique is in principle able to estimate T_eff, log g and [Fe/H] from the AAOmega spectra in the Ca-triplet region, we find in practice that there are too many degrees of freedom at our resolution, which leads to larger errors in the stellar parameters. It is better to estimate T_eff independently. This is consistent with the findings of Kirby, Guhathakurta & Sneden (2008) at their slightly lower resolution of 8000. After correction for the reddening (Schlegel, Finkbeiner & Davis 1998), stellar temperatures are derived from (J − K)_○ using a metallicity-dependent calibration due to Bessell, Castelli & Plez (1998).¹ The procedure is iterative: starting with an initial estimate of [Fe/H] to get an initial estimate of temperature, new estimates of gravity and [Fe/H] are then determined from the χ² technique against a grid of models. Typically four cycles are needed for convergence. The mismatch between isochrones and stellar parameters seen in Fig. 3 for the more metal-rich stars could possibly be due to problems with the stellar models, the adopted reddening or our temperature calibration. It seems unlikely that our abundance scale or our log g scale is the problem: they have both been carefully checked from observations of stars in previously well-observed star clusters covering a wide range of [Fe/H] and log g.

The BaSTI isochrones appear to work well for high-metallicity giants (Cassisi 2010). Although we note that there is some evidence for a spread in age for the metal-rich dwarfs in the bulge (Bensby et al. 2010), which could lead to some displacement towards higher temperatures for the corresponding metal-rich giants, we will assume that the 10 Gyr BaSTI isochrones are appropriate for our bulge stars.

Some authors argue that the Schlegel reddening may be too high by 30–50 per cent for regions of A_V > 0.5 (Arce & Goodman 1999; Cambrésy, Jarrett & Beichman 2005). Although our fields were selected to have relatively low reddening, the regions of highest reddening in our survey which are closest to the Galactic plane also correspond to the locations of the most metal-rich stars. On the other hand, (1) our stars all lie at Galactic latitude |b| > 4°, in the inner Galaxy. These stars are at least 300 pc above the Galactic plane, so the Schlegel reddening should be appropriate; (2) when selecting fields, we found that the Schlegel reddening was in excellent agreement with the reddening inferred from the morphology of the NIR colour–magnitude diagram; and (3) our most highly reddened field at (l, b) = (0°, −5°) has a mean A_V ∼ 2.2. If the Schlegel reddening is systematically too high, then we would expect to see a corresponding shift relative to the isochrones across all metallicity bins in this field. This was not seen. In summary, we believe that systematic overestimation of the reddening is not responsible for the offset seen in Fig. 3 at high metallicities. However, although a systematically high reddening is unlikely to be responsible for the offset seen in Fig. 3 at high metallicities, an overestimation of the reddening for the highest reddening fields only may explain the offset. Fig. 4(a) shows the position of the stars on the isochrone for the most metal-rich stars only ([Fe/H] > 0.25), where the temperature has been calculated for a reddening of 0.9 × E(J − K) only for E(J − K) > 0.3, where E(J − K) denotes the Schlegel reddening. The metal-rich stars now lie closer to the isochrone and the other metallicity bins are not affected. The effect on the derived MDF of using these lower temperatures is that the number of metal rich stars with [Fe/H] > 0 would decrease by about 15 per cent and the mean [Fe/H] of the MDF would decrease by 0.06 dex.

We performed our analysis prior to the availability of the Gonzalez et al. (2011) reddening map, which provides only partial coverage of our fields. The Gonzalez et al. (2011) reddening values in our fields are lower than the Schlegel values. For our field, (l, b) = (0°, −5°), the Gonzalez maps are about E(J − K) = 0.1 on average lower in estimated extinction compared to the Schlegel maps. This reduces the determined temperature by about 300 K (which corresponds to a lower mean [Fe/H] by about 0.2 dex). For the highest metallicity bin (Fe/H > 0), this lower photometric temperature gives a better fit to the stars on the isochrone, but at all other metallicities (Fe/H < 0), the stars are shifted too low in temperature and off the isochrone.

The temperature calibration is the other possible source of the offset. Our adopted temperatures have been used for the calibration of our parameters against a number of well-measured stars and clusters (see ARGOS II), and this temperature calibration is able to return the log g and [Fe/H] of our comparison stars within |$\sigma _{\text{log}\,g}$| = 0.3 and |$\sigma _{\text{[Fe/H]}} = 0.09$|⁠, respectively. The temperature relationships selected were for α-enhanced models for [Fe/H] ≤ −0.25 and for solar scaled models for [Fe/H] > −0.25. Over the T_eff range 4000–6000 K, a change of 50 K in T_eff corresponds to a change of 0.05 in [Fe/H] and 0.15 in log g. A systematic temperature error as large as the offset seen in Fig. 3 seems to be excluded (unless there is a systematic problem with the grid of model spectra such that the model spectra give the right [Fe/H] and log g values for stars whose temperatures have been overestimated from the photometry). We note however that over our temperature range, the Alonso, Arribas & Martínez-Roger et al. (1999) temperatures derived from (J − K) are typically up to 180 K cooler than the scale which we have adopted; the difference is largest for the more metal-poor stars. Fig. 4(b) shows the improved fit to the isochrone by adopting the α-enhanced temperature relationship from Bessell et al. (1998) also for the metal-rich stars ([Fe/H] > −0.25). The cooler temperatures move the metal-rich stars closer to the isochrone. The effect on the derived MDF using these lower temperatures is that the total number of stars with [Fe/H] > 0 decreases by 15 per cent and the mean [Fe/H] of the MDF decreases by 0.03 dex. In summary, although the adopted temperatures give accurate [Fe/H] and log g values, the temperature scale remains a possible reason for the offsets relative to the isochrones. Note that for the T_eff–(J − K) relation, the J − K changes relatively slowly with T_eff, so it is at more risk of errors than an alternative colour index such as V − K. However, our stars are sourced from the 2MASS catalogue for which J − K is available and this colour relation is also less sensitive to differential interstellar reddening which may be a problem for a few of our b = −5° fields.

Following this discussion, we formulate the following procedure for estimating stellar distances (see also Paper II). For stars that are located near the clump (log g = 1.8–3.2 and T_eff = 4500–5300 K), we adopt an absolute magnitude of M_K = −1.61 (Alves 2000) which is appropriate for stars with [Fe/H] > −1.0 (95 per cent of our sample). See Appendix A for further discussion of the absolute magnitude of the red clump adopted for distance calculations. As can be seen from Fig. 2, the majority of stars are taken to be clump stars. For stars outside this range of T_eff and log g, we derive M_K values using the BaSTI isochrone of age 10 Gyr and appropriate [Fe/H], adjusting only the temperature of the star to place it on the closest branch of the isochrone.

For the red clump stars, there are three main sources of distance uncertainty: (i) the Schlegel reddening uncertainty, (ii) the 2MASS uncertainty in the K magnitude and (iii) the red clump intrinsic magnitude dispersion. The reddening uncertainty of the Schlegel maps is estimated to be 16 per cent (Schelgel et al. 1998). The mean reddening in our fields is E(B − V) = 0.2. A 15 per cent reddening error would typically contribute |$\sigma _{K_{S}}$| = 0.02 to the K magnitude. In our stellar selection, we selected only those stars with 2MASS photometric errors ≤0.06 in each band. The mean error in each magnitude for the stars in our sample is much smaller, about 0.03 mag.

The red clump magnitude dispersion estimated by Alves (2000) is |$\sigma _{K_{A}}$| = 0.22. Although this value includes the intrinsic error of the parallax measurement which has |$\sigma _{M_{K_\pi }} = 0.1$| (Alves 2000), we adopt the intrinsic dispersion of the clump to be σ_K ≈ 0.22 to include the small contribution from reddening and magnitude errors described above. As discussed in Paper II, in the (log g, T_eff) plane, the red clump giants lie on a background of first ascent giants. The first giants with log g values close to the red clump have M_K values close to that of the clump. However, because our measuring errors for log g are relatively large (0.30 dex), we have taken clump stars to be those with log g between 1.8 and 3.2, and our sample of clump giants therefore includes some brighter and fainter stars from the first giant branch (see Paper II, Fig. 7). For an unbiased population, most of these first giant branch stars are fainter than the clump stars. Our sample is not unbiased, because the 2MASS catalogue is becoming incomplete for our fainter stars, and the numbers of first giant branch stars that we have included as clump stars are therefore reduced relative to the unbiased sample.

For stars not on the clump, the distances are determined from the temperature and gravity via the BaSTI (Cassisi et al. 2006) isochrones. To determine its absolute magnitude, the star is shifted across in temperature to fit the nearest branch of isochrone. The error on log g (±0.30) dominates the error of the absolute magnitude. For the 10 Gyr BaSTI isochrones, the absolute K magnitude increases almost linearly with log g, with dK/d(log g) = 2.35. Therefore, |$\sigma _{M_K} = 2.35 \times 0.30 = 0.70$| for stars not on the clump. This corresponds to a distance error of about 38 per cent.

By using the red clump membership for as many of our potential clump stars as possible, across 1.8 < log g < 3.2 we are able to take advantage of the precision of the red clump absolute magnitude in calculating distances. Selecting stars over this interval has the advantage of including 99 per cent of our red clump stars given our error of 0.3 dex in log g. The compromise in doing this is the contamination of non-clump stars that will be intrinsically brighter or fainter which will be brought into the inner region. For example, stars with a true log g at the upper limit of this selection will in reality lie closer to us than the bulge region which we define in our analysis to be R_G < 3.5 kpc. These stars will be moved further in distance from the Sun into the bulge. Similarly, stars with a real log g at the low end of this selection will in reality be at larger distances than the bulge region. This kind of contamination by first ascent giants will be a problem for any distance estimates using the red clump magnitude (e.g. Rattenbury et al. 2007), unless gravities are accurately known (Valentini & Munari 2010). Fig. 5 shows the distribution of the stellar distance to the Sun as a function of log g for stars taken with distances from the isochrone only (at bottom) and using the red clump magnitude across 1.8 < log g < 3.2 (at top) for our minor axis fields. These figures demonstrate how taking a range in log g as assigned clump stars to allow for our observational errors moves a fraction of closer and farther stars into the bulge region. The boxy structure seen in Fig. 5 across xy reflects the limits across our log g range that we have assumed for the clump and the K-magnitude range selected for our observations of 11.5–14. The extent of the spread introduced in taking a range for the clump along the y-axis is set by this magnitude range. We can reduce only our interval in log g, along the x-axis. We could have used a narrower interval in log g of say 2.0–3.0 but then we would have lost about 13 per cent of the clump stars. Fig. 6 evaluates the effect of our selection on the MDF for stars within R_G < 3.5 kpc on the minor axis using this interval in log g for red clump membership. The MDF at left is for the population where all distances have been taken directly from the isochrone. There are about 10 000 stars in the inner region using isochrone distances only. The centre panel shows the MDF generated for stars calculated with isochrone+red clump distances and there are about 14 150 stars in the inner region using isochrone+red clump distances. The MDFs are only weakly affected. The far-right panel shows the fractional difference of these distributions of MDF (isochrone+red clump)–MDF (isochrone only). The selection difference is small in the inner 3.5 kpc, and taking selections of smaller radii increases the relative number of metal-rich stars and decreases the number of metal-poor stars. For a selection cut of R_G < 1.5 kpc, there are about 10 per cent fewer stars in total with −0.2 > [Fe/H] > −1.0 and 10 per cent more stars with −0.2 > [Fe/H] > +0.5.

To summarize, we take all of our potential clump stars with log g between 1.8 and 3.2 to have M_K = −1.61 ± 0.22 (about 900 pc at the distance of the Galactic Centre), recognizing that for some of the underlying first ascent giants the uncertainty in M_K will be somewhat larger. For the non-clump giants, the distance error is 0.70 mag, dominated by the uncertainty in the log g estimates. In our data analysis, we verified that similar results were obtained for the smaller sample selection obtained using isochrone-only distances. In Section 3, we evaluate Gaussian component fits to our MDFs. The relative fractions of the components are only slightly changed if we use isochrone-only distances.

To evaluate the properties of the stars which lie in the inner Galaxy, we take stars with cylindrical Galactocentric |R_G| ≤ 3.5 kpc to be part of the inner Galaxy. This sample includes about 14 150 of the original 25 500 stars. The metallicity distribution of these inner stars is summarized in Table 1. They are a mixture of all of the populations present in the inner Galaxy: the metal-poor halo, the thick disc, the thin disc around the bulge and the bulge itself.

2.2 The line-of-sight distribution of giants in the inner Galaxy

Our sample includes giants all along the lines of sight with distances up to ∼55 kpc from the Sun on the far side of the bulge. The observed distribution of stars over distance is affected by our selection criteria, observational strategy and the stellar metallicities. Fig. 7 shows the distribution of radii R_G for our entire sample of 25 500 stars. This figure does not reflect the underlying density distribution. Negative radii are for stars on the far side of the bulge, and the Sun is located at R_G = 8 kpc. The observed stellar density distribution is shown as the smooth curves (generalized histograms with a kernel of σ = 1.5 kpc) for (i) all stars and (ii) stars in intervals of [Fe/H] to compare the relative density distribution of different metallicity groups to the total distribution of stars. In total, 60 per cent more stars are observed per kpc along the line of sight on the near side of the bulge than on the far side. However, the line-of-sight distribution is metallicity dependent. The majority of observed giants are located about 4–10 kpc from the Sun, with the peak about 3 kpc in front of the bulge. Robin et al. (2012) estimate a peak in the thin disc number density at ∼2.3 kpc from the centre of the Galaxy which may contribute to the larger number of metal-rich stars seen at low Galactic latitudes on the near side of the bulge. The more metal-poor stars are less biased towards the near side of the Galactic Centre, with only 15 per cent fewer stars observed in this range on the far side of the bulge. We preferentially sample the more metal-poor stars at larger distances because they are bluer, brighter and located at higher galactic latitudes.

Fig. 7 shows a minimum or inflection in the stellar number density around R_G = 0. This is due to the location of our fields; stars with |R_G| < 1.0 kpc come only from our three fields at longitude l = 0° and at latitudes b between 4° and 11°. This central feature is much less prominent for the more metal-rich stars, indicating that the more metal-rich stars are more concentrated to the inner regions of the bulge. The upper panel of Fig. 7 shows that most stars with |R_G| < 0.75 kpc have [Fe/H] > −1.0: only two stars from our sample in this region have [Fe/H] < −1.0. We will return to these properties of the inner region when we examine the MDF as a function of radius R_G (Figs 12 and 13).

Although the lower limit of metallicity for our stars depends on radius (higher for |R_G| < 0.75 kpc), the upper limit of the [Fe/H] distribution is similar at about [Fe/H] = 0.4 for all radii up to about |R_G| = 5 kpc, with only a few stars above this level. Table 2 shows the distribution of the number of stars (giants) with R_G for stars seen out in the halo on the far side of the bulge.

The most metal-poor stars in our sample, with [Fe/H] < −2.0, are located in the inner Galaxy but not in the innermost region. Of the 27 stars with [Fe/H] < −2.0, 18 lie between |R_G| =0.9 and 4.5 kpc. In the inner region, our most metal-poor star has [Fe/H] = −2.6 and lies at R_G = −3.8 kpc on the far side of the bulge. Our most metal-rich star with [Fe/H] = 0.68 is at R_G = 2.0 kpc on the near side of the bulge. Although we find a few stars at metallicities [Fe/H] ∼ 0.6, which are among the highest [Fe/H] values for bulge stars reported in the literature, we note that our χ² abundance scale is externally calibrated only for 0.1 > [Fe/H] > −2.1. Outside the range, the abundance scale comes directly from the grid of model spectra.

We have not included the solar neighbourhood dwarfs in Fig. 7. In total, our sample has about 2000 dwarfs located between about 1 and 3 kpc from the Sun. They have a median [Fe/H] = −0.1. These stars are not relevant to the present discussion of the inner Galaxy and are not included in the discussion here.

2.3 The [Fe/H] distribution in the inner Galaxy

From the mean of the abundance distributions in our three latitude zones along the minor axis, we confirm the metallicity gradient seen in the bulge by Zoccali et al. (2008b). Our MDFs in the three minor axis fields include a total of about 2000 stars with |R_G| ≤ 3.5 kpc (see Fig. 8), which is about 70 per cent of the whole sample of stars observed on the minor axis. Fig. 9 shows the MDF for all of our survey stars with |R_G| ≤ 3.5 kpc. This distribution shows the relative proportion of stars in each [Fe/H] bin; it includes about 14 150 stars with contributions from fields across longitudes between ±20° and all four latitudes of b = −5°, −7| ${.\!\!\!\!\!\!^{\circ}}$|5, −10° and b = +10°.

Our inferred MDFs shown in Fig. 8 are qualitatively similar to those obtained by Zoccali et al. (2008b) from their survey of about 800 K – giants observed at high spectral resolving power (R ≈ 45 000) in their minor axis fields at b = −4°, −6° and −12°. However, our selection includes a larger fraction of metal-poor stars due to our significantly bluer colour cut in our stellar selection compared to that of Zoccali et al. (2008b). Quantitatively, our mean [Fe/H] values are about 0.15 dex lower on average, with our mean [Fe/H] decreasing from −0.16 at b = −5° to [Fe/H] = −0.5 at b = −10°. We find a somewhat lower vertical gradient of the mean abundance along the minor axis: −0.45 dex kpc⁻¹ compared to their value of −0.6 dex kpc⁻¹.

3 THE COMPONENTS IDENTIFIED IN THE [Fe/H] DISTRIBUTION

We now examine the structure of the MDFs in our three latitude zones. To eliminate the effect of binning we examine the generalized histograms for our data in our fields, where each star is represented by a Gaussian distribution centred on its [Fe/H] value and with a standard deviation of 0.13 (we adopted a kernel that is marginally wider than our estimated measuring error). To increase the number of stars, we have binned our data across all of the bulge fields, out to l = ±15°, for each latitude b = −5°, −7| ${.\!\!\!\!\!\!^{\circ}}$|5 and −10°. The fields in the three latitudes which are grouped together are shown in Fig. 10. For a cut in radius of |R_G | ≤ 3.5 kpc, there are 4200 stars in the lowest latitude fields (b = −5°), 2000 stars in the intermediate-latitude fields (b = −7| ${.\!\!\!\!\!\!^{\circ}}$|5) and 4000 stars in the highest latitude fields (b = −10°). It was immediately obvious that the MDFs showed some structure which was qualitatively similar over our three latitude intervals.

Figs 11(a), (b) and (c) show the generalized histograms at these three latitudes for stars within the bulge region |R_G| ≤ 3.5 kpc. The histograms show some structure and we can decompose them into a few Gaussian components. For all of the three latitudes, we can fit almost the same four Gaussian components, which we denote by A, B, C, D² plus an additional more metal-poor component E in the higher latitude fields (E is almost negligible at b = −5° and comprises only a very few stars at all latitudes). Figs 12(a), (b) and (c) show conventional histograms of the metallicity distribution in the three latitude zones. The individual MDF components can be identified. The smooth curves are Gaussians with the same mean [Fe/H] as for the components in the generalized histograms of Fig. 11, but with standard deviations |$\sigma _{\text{real}}$| reduced to account for the broadening effect of the generalized histogram kernel (σ_k = 0.13): the |$\sigma _{\text{real}}$| values for the Gaussians shown in Fig. 12 are given by |$\sigma ^2_{\rm real} = \sigma ^2_{\text{GH}} - \sigma ^2_k$|⁠, where |$\sigma _{\text{GH}}$| are the standard deviations derived in the multi-Gaussian fit to the generalized histograms.

Table 3 summarizes the parameters for the MDF fits, corrected for the standard deviation of the generalized histogram kernel. We find slightly decreasing mean abundances of components A, B and C with increasing height below the plane |z|.³

Gaussian fits to the components have been adopted here for simplicity. It may be more correct to use skewed distributions for the components (e.g. Gilmore, Wyse & Jones 1995). This may remove the need for the extra Gaussian components at the low-metallicity tail of the distribution, for which there are only a few stars. Components D and particularly E in the metal-poor tail are fitted here only to demonstrate the changing fraction of metal-poor stars. Where stars with Fe/H < −1.5 are present, component E is included to highlight their relative fraction, although this component is generally not statistically significant. It is instructive however to see how this metal-poor contribution changes with latitude and longitude. Table 3 demonstrates how the fraction of metal-poor stars increases with latitude from an almost negligible value in our lowest latitude zone.

A Monte Carlo analysis was performed to determine the error in the MDF component parameters due to the fitting and sampling. We took the observed sample of [Fe/H] values and generated a new value of each [Fe/H], drawn from a Gaussian distribution around the mean of the data point and its standard deviation (we adopted σ = 0.13, fractionally larger than the measuring error). Revised Gaussian fits were then made to the generalized histogram generated from this new distribution and the process was repeated 10³ times to determine the variance in the mean, amplitude and standard deviation of the five components. Table 3 summarizes the outcome: for the mean [Fe/H] of the three main components A, B, C shown in Figs 11(a)–(c), the errors are ≈0.01. For the amplitudes they are typically 0.01, and for the standard deviations of the Gaussian components, they are about 0.01–0.05. The corresponding Gaussian components along the minor axis fields only are shown in Table 4.

3.1 Justification of components

The maximum likelihood method was employed using the Expectation–Maximization (EM) algorithm of Meng & Rubin (1993) for multi-Gaussian fits within the statistical package r,⁴ to determine the log likelihood for N = 2, 3, 4, 5, 6 components for our MDF distributions given our probability distribution function.

The BIC method of assessing goodness of fit, which uses the maximum likelihood (and associated Akaike information criterion) is simplified metrics from a fully Bayesian analysis (Heavens 2010). This method finds only a maximum and not average likelihood, however. The complete Bayesian approach calculates the likelihood over the full parameter space, thus assessing the predictive power of the model, given a prior. However, for our purposes, the information theory approach is sufficient to find the best model. This is because although we aim to demonstrate that our fits are valid, our decomposition is also motivated by the peaks we can see in the smoothed distribution, the behaviour we find for the red clump stars along the minor axis (Ness et al. 2012) and the clear dependence of kinematics on [Fe/H] (Ness et al., in preparation).

The BIC test returns a minimum (optimal) value for a four-component fit (A, B, C, D) at b = −5° and −7| ${.\!\!\!\!\!\!^{\circ}}$|5 and a minimum for a five-component fit (A, B, C, D, E) at b = −10°. In Figs 11, 13 and 14 and Tables 3 and 4, the component E at the lowest [Fe/H] group is included to demonstrate the changing number of the most metal-poor stars seen in the sample as a function of latitude, even though this is not a statistically significant fit. The maximum likelihood fitting routine returns the parameters of the χ² minimization fitting of Table 3 (except at the metal-poor tail) for the MDFs at b = −7| ${.\!\!\!\!\!\!^{\circ}}$|5 and −10°. For these mixture distributions, component D is returned with weighting ≈ 0.01 and with a wider standard deviation than for the χ² fit (σ ≈ 0.5). This discrepancy does not affect the other Gaussians fits nor does it alter the conclusions from the decomposition. The agreement of the two approaches indicates that the majority of stars are in components A, B and C and these components are securely fitted. The metal-poor components D and E involve a relatively few stars and their contributions for the different fields can be seen in Table 3. They are not the main target of our discussion here.

3.2 Consequence of mixture model

The component fitting described above identifies the three major Gaussian components A, B and C which are directly visible in the MDFs, plus the two minor components D and E. We are motivated here not only by these fits but also by the split clump (see Ness et al. 2012) and kinematics which are clearly different for our identified components (see Ness et al., in preparation). Independent of our MDF component fits, when we examine the kinematics and distribution of the stars partitioned into [Fe/H] intervals of 0.5 dex from +0.5 to −1.5 (these intervals correspond approximately to our components A–D), these subgroups show different kinematics and different split clump characteristics.

If we were to adopt two components for our data rather than our best-fitting representation of four or five components, the optimal two components comprise a narrow metal-rich distribution and a wider metal-poor distribution, with (μ₁, σ₁) = (0.15, 0.1) and (μ₂, σ₂) = (−0.3, 0.35). This is similar to the components identified by Hill et al. (2011) and by Babusiaux et al. (2010) in their b = −4° field. Both argue for a two-component bulge population, comprising a metal-rich population and an old spheroid. The number of components one considers clearly plays an important role in the interpretation.

4 PROPERTIES OF THE COMPONENTS

The decomposition of the MDF into the A, B, C (D and E) components allows us to quantify the varying contributions of the metallicity components over l and b and also over distance along the line of sight.

Fig. 11 shows how the relative contributions of the components change with latitude. At b = −5°, the more metal-rich components A and B dominate, with a smaller contribution from component C and only a few stars from the metal-poor components D and E. At b = −7| ${.\!\!\!\!\!\!^{\circ}}$|5, components B and C are the strongest; component A is clearly weaker than at b = −5°. At b = −10°, components B and C still dominate and component A is now quite weak. The metal-poor components increase in significance with increasing | b |, and more metal-poor stars ([Fe/H] < −1) appear at higher latitudes.

In the discussion in Section 7 about the nature of the metallicity components seen in the inner Galaxy, we show MDFs for intervals of R_G which demonstrate the changing contribution fraction of the different components with distance from the Galactic Centre. Fig. 13 shows the MDF components in various distance intervals as indicated below the individual panels for b = −5°, and Fig. 14 is similar for b = −10°. The data for each latitude zone are combined as before for l = ±15°.

Figs 13 and 14 show that the components A–E are present at both latitudes out to |R_G| ≈ 5 kpc. The widths of the individual components are approximately unchanged, but their contribution fraction changes with latitude and R_G. Examined separately, the near and far sides along the line of sight show comparable MDFs out to |R_G| < 4.5 kpc. There are fewer metal-rich stars observed on the near side however and a radial change in the morphology of the MDF occurs for |R_G| > 4.5 kpc: at these larger radii, components A and B merge (see for example components A and B in panels e and f of Fig. 13). Note that this conclusion depends on the adopted metallicity gradient for component A as discussed below. The distance cuts show that within R_G < 4.5 kpc the two components A and B are necessary to fit to the stars with [Fe/H] > −0.5. However, outside the bulge (for a radius of R_G > 8 kpc at b = −5° and R_G > 4.5 kpc at b = −10°), a single component is sufficient to fit to the stars with [Fe/H] > −0.5 (see Figs 13i and 14i).

At b = −5°, component B dominates between 0.75 and 3.0 kpc but inside of this region the most metal-rich component A is marginally strongest. Out to 4.5 kpc, component C is no more than about 50 per cent of that of A. At b = −10°, component B is the strongest out to |R_G| = 6 kpc and component C is a larger fraction of component B out to 4.5 kpc compared with the lower latitude. Within R_G < 4.5 kpc, component A is relatively weak, at no more than about 30 per cent of the strength of C.

At b = −10°, component A which appears as a small fraction in the inner region is no longer present on the far side of the bulge (see Figs 14g and h). At both latitudes, component B appears in the largest fraction on the far side of the bulge and component A is the largest fraction on the near side of the bulge outside of |R_G| > 4.5 kpc.

Components D and E include only a small number of stars. At b = −5°, their fraction increases with Galactic radius. The fraction of very metal-poor stars (in D and E) is higher at b = −10°, particularly near the Galactic Centre (|R_G| < 0.75 kpc).

4.1 Comparison with disc fields

One of our fields lies out in the disc at (l, b) = (−31°, −5°), for which the minimum value of |R_G| is about 4 kpc. The same components A, B and C are seen in this disc field (Fig. 15). From Figs 13 and 14, these three components are also seen in the radial intervals extending out into the disc at smaller | l |. The presence of the same three components A–C in the bulge and the inner disc makes sense. In the inner disc, these components are present in their undisturbed form, as disc stars, while in the bulge fields the stars of these components have been incorporated dynamically into the bulge.

With increasing distance from the Galactic Centre, we note again the change in component A relative to components B and C. For |R_G| > 4.5 kpc, components A and B merge, and at b = −10°, four components (B, C, D and E) are sufficient to fit the data, assuming a radial gradient in the mean abundance for A. At b = −5°, where components A and B merge outside of the inner bulge region, component A is dominant close to the Sun. On the far side at R_G > 4.5 kpc, component B is dominant and components A and C are roughly equal out to R_G = 7.0 kpc. At b = −10°, component A is similarly the most significant near the Sun and at this higher latitude on the far side of the bulge, for R_G > 4.5 kpc, component C is the strongest. Table 5 shows the relative heights below the plane for each radial cut. We note that stars with a distance of about 6 kpc on the far side of the bulge at b = −10° are located at z-heights of 2.4 kpc below the plane. At b = −10° and |R_G| > 7 kpc, the stars are on the far side of the bulge at mean z-heights > 2.6 kpc. Therefore, it is not surprising that the component C (with its thick-disc-like metallicity range) should dominate on the far side of the bulge.

4.2 The metallicity gradient with distance |z| from the plane

Along the minor axis, the mean [Fe/H] for some of the individual components of the MDF shows a weak decrease with latitude. From b = −5° to −10°, the total decrease in the mean abundance of components A and B is −0.06 ± 0.04 and −0.06 ± 0.03 dex, respectively. The mean abundance of component C shows no significant decrease with latitude on the minor axis, although the longitude-integrated MDF does show a small decrease in the mean of component C with latitude (see Table 3) of −0.03 ± 0.02 dex. The numbers of stars with [Fe/H]<−1 are small (10, 29 and 55 stars at b = −5°, −7| ${.\!\!\!\!\!\!^{\circ}}$|5 and −10°, respectively, along the minor axis) and we cannot evaluate the abundance trends with height for the metal-poor components D and E.

A small decrease in the mean abundances of components A and B for the integrated MDFs (l = ±15°) is also measured for these components. This small change in the mean abundances with latitude, and the large change of the contribution fractions of the components, indicates that the overall observed abundance gradient measured by Zoccali et al. (2008b) from b = −4° to −12° (of −0.6 dex/kpc between −4° to −6°) comes mainly from the changing proportions of the MDF components with latitude, as shown in Fig. 11, rather than from gradients in the individual components.

From b = −5° to −10°, the gradients for A and B with height below the plane calculated from the minor axis MDFs are −0.08 ± 0.05 dex kpc⁻¹ and −0.08±0.04 dex kpc⁻¹, respectively. No gradient in C is seen along the minor axis but the integrated MDFs show a small gradient of −0.04 ± 0.02 dex kpc⁻¹.

4.3 The radial metallicity gradients

The mean abundances for components A, B and C as a function of R_G are given in Table 5 with their corresponding errors. A weak radial gradient (≈ −0.02 dex kpc⁻¹) is evident in the mean abundance of the most metal-rich component A. We note that this gradient may come partly from a radial gradient and partly from a z-gradient in these zones of constant latitude. This weak abundance gradient is shallower than the gradient seen in the cepheids (Luck & Lambert 2011) and young (age ≤ 1 Gyr) open clusters (Yong, Carney & Friel 2012). However, given an additional negative gradient with height from the plane, this measured gradient will be smaller than the true value. This is because with increasing distance along the line of sight, the stars sampled are further from the plane and will be relatively more metal poor.

At R_G > 8 kpc at low latitudes and R_G > 4.5 kpc at high latitudes, the most metal-rich fraction of stars is no longer present and the MDF can be described by two main components: the more metal-rich component is component B. This component represents the most metal-rich stars that we sample at large heights from the plane; we note that at b = −5° the |z|-height at 6 kpc on the far side of the bulge is 1.2 kpc. From the inner region, component B shows almost no abundance gradient with |R_G|. The abundance of component B is higher on the far side of the bulge where components A and B have merged (b = −5°) or component A is no longer present (b = −10°).

The mean abundance of component C remains relatively flat at both latitudes out to about |R_G| = 4.5 kpc. Outside this radius, the errors are too large to evaluate any gradient. At R_G > 8 kpc, the |z|-height above the plane is about 3 kpc for b = −10° and there may be increasing contamination of component C by the metal-poor stars which belong to the overlapping lower metallicity component D.

We have used the above gradient in component A to argue in the previous section that components A and B merge at a radius of |R_G| ≈ 4.5 kpc. Within R_G < 4.5 kpc at both latitudes, the MDF shows three significant peaks: two at [Fe/H] > −0.5 and one at −1.0 > [Fe/H] > −0.5. Outside of R_G > 4.5 kpc, there is only one clear peak in the MDF present at [Fe/H] > −0.5 and a second significant peak at −1.0 > [Fe/H] > −0.5. At all radii, there is the small contribution of stars with [Fe/H] < −1.0 represented by components D and E.

5 α-ENHANCEMENT

The α-enhancement reflects the rate of chemical evolution. Stars are α-enhanced if the chemical enrichment in their formation environment has been sufficiently rapid that there was not time for the Fe-enriching SN Ia to contribute fully to their chemical evolution. In the following, [α/Fe] denotes the mean of [Si/Fe], [Mg/Fe] and [Ti/Fe] (we used five Si lines, three Mg lines and nine Ti lines); we were unable to make an accurate estimate of [Ca/Fe]. Fig. 16 compares the radial distribution of α-enhancement for stars in 0.5 dex bins of [Fe/H] from −2.0 to 0.6, extending in radius from inside the bulge to the outer disc and halo at R_G = 30 kpc. For this figure, we have not associated the [Fe/H] intervals with the MDF components, because of the abundance gradient seen in the most metal-rich component A.

5.1 α-enhancement trends with Galactic radius

Stars with [Fe/H] > −0.5, and particularly the most metal-rich stars, show a decrease in [α/Fe] from slightly enhanced values in the inner Galaxy. Stars with [Fe/H] > 0 show [α/Fe] declining to negative values beyond |R_G| > 20 kpc (Fig. 16). This decrease in the [α/Fe] ratio with radius indicates that the chemical evolution in the inner Galaxy proceeded more rapidly than in the disc at larger radii where the [α/Fe] ratio is near solar for the most metal-rich stars. Note that stars at these larger values of |R_G| > 7 kpc are mostly on the far side of the bulge, at z-heights >1.4 kpc.

For stars with −2.0 < [Fe/H] < −0.5, the α-enhancement is fairly flat with galactic radii out to 30 kpc. The most metal-poor fraction of stars shows a slight increase in enhancement, although these stars show a larger scatter.

Fig. 17 shows the radial distribution of mean [α/Fe] for the most metal-poor stars with [Fe/H] < −2.0. (The most metal-poor star in the sample is at [Fe/H] = −2.9; however, it is difficult to measure the gravity from the weak Ca-triplet lines and we have no distance estimate for this star.) The distribution of [α/Fe] for these metal-poor stars shows a steep decrease from the inner region. In this subsample of metal-poor stars with |R_G| < 8 kpc, 75 per cent lie within |R_G| < 4.0 kpc.

5.2 α-enhancement and components

We can analyse the α-enhancement in the context of the components of the MDF. Fig. 18 shows contour plots of the stellar density for stars with R_G < 3.5 kpc in the [Fe/H]–[α/Fe] plane, within longitude ±15°, with the generalized MDFs for the components A, B and C superimposed. At b = −5°, there are two clear peaks in the density distribution, corresponding to components A and B. At b = −10°, component A has attenuated relative to component B and the contour plots show only a single density maximum associated with component B. The most metal-rich component A has a mean [α/Fe] of about 0.1 at b = −10° and 0.16 at b = −5° [in agreement with the weakly enhanced [α/Fe] at solar [Fe/H] found by high-resolution studies of bulge stars in Baade's window by Alves-Brito et al. (2010)]. All of the lower metallicity components in the bulge are more α-enhanced. The α-enhancement increases as expected with decreasing [Fe/H], from about 0.1 at solar metallicity to about 0.6 at [Fe/H] = −2.0. The dispersion in [α/Fe] increases as [Fe/H] decreases. At latitude −10°, the contribution of the more metal-rich component of stars is attenuated, and we see a lower dispersion in the [α/Fe] for [Fe/H] > 0 relative to the lower latitude strip. The χ² error on the [α/Fe] measurement is σ = 0.1 and this is consistent with the spread seen in Fig. 18.

Fig. 19 compares the [Fe/H]–[α/Fe] relation for stars in three radial zones: inside the bulge at |R_G| < 1.2 kpc and |R_G| < 3.5 kpc, and outside the bulge at 5 < |R_G| < 10 kpc. For the innermost interval inside the bulge region, the α-enhancement increases linearly with decreasing [Fe/H], from about 0.1 for stars above solar metallicity to about 0.6 at [Fe/H] = −2.0. The dispersion in [α/Fe] increases more slowly, from about 0.1 at the most metal-rich end to about 0.15 at [Fe/H] = −2.0. This dispersion at the more metal-poor end is higher for stars within |R_G| < 3.5 kpc, at σ ≈ 0.35; other dispersion values in this radial interval follow those of the inner region R_G < 1 kpc.

Comparing the α-enrichment of stars in the disc and the bulge, we see that the stars with 5 < |R_G| < 10 kpc and [Fe/H] > −0.5 have similar [α/Fe] values to those inside the bulge. The [α/Fe] values are about 15 per cent lower for stars with [Fe/H] < −1.0 and about 10 per cent higher for stars with [Fe/H] > −1.0 relative to those with |R_G| < 3.5 kpc. The dispersion is higher than for |R_G| < 3.5 kpc by about 10 per cent at the more metal-rich end and 30 per cent at the more metal-poor end.

6 INTERPRETING THE COMPONENTS

In this section, we discuss two significant features of the MDF components: (1) their individual vertical abundance gradients and the overall vertical abundance gradients associated with the changing fraction of the major components with height above the plane, and (2) the contribution of the components to the split red clump of the bulge.

6.1 The vertical metallicity gradient of the bulge

In the bulge region, components A, B and C individually show weak vertical abundance gradients (Section 4.2) from −0.08 ± −0.05 to −0.04 ± 0.02 dex kpc⁻¹ over the latitude interval b = −5° to −10°. Although our overall vertical metallicity gradient (−0.45 dex kpc⁻¹) is comparable to that found by Zoccali et al. (2008b), we interpreted this overall gradient as coming from changing fractions of components A–D with latitude via the action of the instability process on the pre-existing early disc.

It is not clear how to interpret these small vertical abundance gradients in the individual components. Are they the remnants of abundance gradients that were already in place in the thin and thick discs before the instability event? The small gradients which we observe compare well with the fairly flat vertical abundance distribution found by Bekki & Tsujimoto (2011), from their simulation of an initial disc system which buckles into a boxy/peanut bulge/bar via the disc/bar instability. They model both the instability of an initial single thin disc and also a dual (thin+thick) disc population. These discs initially have vertical abundance gradients. Due to the vertical mixing of the stellar populations by the bar, they find that the bulge is unlikely to have a steep vertical metallicity gradient. Their model does not include chemical evolution so is appropriate for a system in which the chemical evolution took place in the early disc, before the instability events and the formation of the bar/bulge.

6.2 The split red clump of the Galactic bulge

There is a striking relationship between the spatial distribution of stars in the bulge region and [Fe/H] which was seen in our analysis of the split clump (Ness et al. 2012). The split clump is a consequence of the peanut structure of the bulge, so the presence of a split clump indicates membership of the peanut bulge.

Fig. 20 shows the magnitude distribution of the red clump stars along the minor axis for our three latitudes, selected with K₀ = 12.38–13.48 and log g = 1.9–3.1. The clump stars are split into three [Fe/H] bins: (i) [Fe/H] > 0, (ii) 0 > [Fe/H] > −0.5 and (iii) [Fe/H] < −0.5. In Ness et al. (2012) we compared the spatial distribution and kinematics of the clump stars with predictions from an evolutionary N-body model of a bulge that grew from a disc via bar-related instabilities. The density distribution of the peanut-shaped model is depressed near its minor axis, thus producing a bimodal distribution of stars along the line of sight through the bulge near its minor axis, very much as seen in our observations for stars [Fe/H] > −0.5. We found that only stars with [Fe/H] > −0.5 show the split structure at latitudes b < −5° and concluded that the stars of the boxy/peanut bulge have [Fe/H] > −0.5. Stars with [Fe/H] < −0.5 do not show the split clump at any latitude. This indicates that the initial disc from which the bulge formed had a very few stars with [Fe/H] < −0.5, or that stars with [Fe/H] < −0.5 were not in a region of phase space which is mapped into the boxy structure by the disc instability.

From Fig. 20, stars with [Fe/H] > −0.5 show the split at b < −5°. Stars in our component A ([Fe/H] > 0) show a deeper split than stars in our component B (0 > [Fe/H] > −0.5) at the intermediate latitude of b = −7| ${.\!\!\!\!\!\!^{\circ}}$|5. At b = −10°, the split is still seen for components A and B but is affected by the asymmetry in the magnitude distribution for the more metal-rich stars: these stars are located preferentially closer to the plane and are therefore found more prominently in the nearer and brighter peak. The proportion of stars with [Fe/H] > 0 decreases from about 40 per cent at b = −5° to about 8 per cent at b = −10°. In contrast, for the component B stars with −0.5 < [Fe/H] < 0, the proportion stays approximately constant from b = −5° to −10°.

The split distribution of the component A and B stars, and the absence of the split for the more metal-poor component C stars, leads us to infer that the peanut bulge structure is defined primarily by the stars of component A ([Fe/H] > 0) and component B (0 > [Fe/H] > −0.5).

7 DISCUSSION

In a survey of the inner Galaxy, we can expect to find stars from the inner regions of all of the major components of the Galaxy. We have identified five metallicity components in total in the stellar MDF of the inner Galaxy and interpreted the MDF as a composite population made up of these components A–E in order of decreasing [Fe/H]. The two more metal-rich components (A and B) define the boxy/peanut structure of the bulge. We make the following interpretation of the metallicity components within the disc instability scenario for bulge formation.

• A: the metal-rich ‘boxy/peanut bulge’ with a mean [Fe/H] ≈ +0.15 that is concentrated towards the plane. In the inner regions, component A has the chemical properties of a younger chemically evolved thin disc (higher [Fe/H] compared to the thin disc in the solar neighbourhood and slightly higher [α/Fe]). At larger Galactic radii, component A appears to be somewhat more metal poor (see Table 5) and we associate it with the thin disc at larger Galactic radii. Its MDF at larger radii is similar to that of the thin disc near the Sun (e.g. Haywood 2008). Although chemically similar to the thin disc, it replicates the line-of-sight kinematics and spatial properties of the thicker boxy/peanut bulge component B in the inner regions of the Galaxy (see Ness et al., in preparation), but at lower heights above the plane (see Section 4). We note again that stars in the metallicity range of component A are prominently involved in the bimodal peanut structure of the bulge seen at higher Galactic latitudes.

• B: the vertically thicker ‘boxy/peanut bulge’ contributes a similar fraction of the bulge stars over the latitude range from b = −5° to −10°. Its mean [Fe/H] value is about −0.25 and it is the dominant component in most of our bulge fields inside R_G < 4.5 kpc (see Figs 13 and 14). In the instability scenario, it would represent the early thin disc when there were a few stars with [Fe/H] < −0.5 dex at the onset of the bar buckling instability. Again, stars in the metallicity range of component B show the split structure of the bulge at higher latitudes. At larger radii, component B is associated with the metal-poor stars of the early thin disc that we sample at large heights from the plane.

• C: we interpret component C as the inner thick disc (mean [Fe/H] about −0.70) which is at most weakly involved in the boxy/peanut structure of the bulge. Component D is tentatively ascribed to the metal-weak thick disc (mean [Fe/H] about −1.2). These components were likely present in the inner regions of the Galaxy (|R_G| < 3.5 kpc) before the instability event, and may have been puffed up by this event. We note again that components C and D do not appear to be part of the split structure of the peanut bulge (Ness et al. 2012). At larger |R_G|, these old components are still present in their relatively undisturbed form.

• E: we interpret this metal-poor component as the stars of the inner Galactic halo.

Components A and B are clearly part of the boxy/peanut structure of the Galactic bulge, as shown by the split red clump discussed in Ness et al. (2012). Models by Athanassoula (2005a) show different levels of boxy/peanut behaviour. For the boxy systems, the isophotes are roughly parallel to the equatorial plane. The isophotes of the peanut systems show a clear dip near R_G = 0. If this dip is relatively deep, the systems give an impression of an X-shaped structure most clear at larger heights above the Galactic plane, particularly after unsharp masking. We have seen in Ness et al. (2012) that component A shows a more prominent central minimum than component B, but it is not clear whether A is intrinsically more X-shaped or if its deeper minimum reflects a similarly shaped structure but with a smaller scale height than the thicker boxy/peanut bulge represented by component B.

7.1 The origin of the components

In the instability scenario, the boxy/peanut bulge formed out of the early thin disc: after about 1 Gyr, the instability of the disc forms a bar which then buckles vertically into the observed boxy/peanut morphology. This scenario gives a guide to interpreting the components of the bulge. We know from the age–metallicity relation in the solar neighbourhood (e.g. Haywood 2008) that the chemical evolution of the solar neighbourhood proceeded rapidly: 10 Gyr ago, the stars of the thin disc in the solar neighbourhood already had a mean [Fe/H] near solar, with a spread from about +0.5 to −0.5. In the context of inside-out formation of the disc, it is likely that the chemical evolution of the inner Galaxy was even more rapid. By the time the bar-forming instability occurred, the chemical evolution of the inner disc was probably already well advanced.

The instability process would affect all of the inner disc stars, and generates a mapping of the stars of the early disc into the boxy/peanut structure. In this way, the bulge preserves a dynamical imprint of the chemical distribution of the thin disc at the time the buckling occurred; the bulge is a chemical snapshot of the MDF of the early thin disc captured in the extended bar. We have seen that two components, A and B, are involved in the boxy/peanut structure of the bulge. Fig. 11 shows that the more metal-rich component A (mean [Fe/H] = +0.15) is more concentrated to the plane than the thicker dominant bulge component B (mean [Fe/H] = −0.25). This more metal-rich component A is chemically consistent with what we would expect for a younger, colder and more metal-rich thin disc in the inner region at the time of the instability event.

Examination of the N-body models shows that the mapping from disc to bulge depends on the location and kinematics of the stars at the time of the instability. Kinematically colder stars in the disc can suffer substantial radial and vertical migration, and can therefore be strongly involved in the boxy/peanut structure. We will discuss this mapping of the disc into the bulge in more detail in a later paper. At this stage, we suggest that component A may simply represent the responsive cold and metal-rich stars of the thin disc in the inner region which were entrained in the potential of the major boxy/peanut component B. The kinematical selectivity of the mapping is likely to explain the presence of the two distinct components A and B in the bulge at the higher latitudes | b | > 5°.

The thick boxy/peanut bulge population B can be interpreted as the product of the bulk of the thin disc (already a few Gyr old at the time of the instability event) which forms a bar and puffs up into a boxy/peanut bulge (Athanassoula 2008). The mean abundance and abundance range of component B are similar to that of the intermediate-age disc stars (3–8 Gyr old) in the solar neighbourhood at the present time (e.g. Haywood 2008). The bulge structure is visible in NIR photometry and star counts out to a radius R_G of about 3.5 kpc. Further from the Galactic Centre, stars of the old thin disc (component B) are still present, but in their relatively undisturbed original thin disc structure.

Component C has a mean abundance like that of the thick disc in the solar neighbourhood. If it is representative of the thick disc in the inner Galaxy before the instability event, it may have been correspondingly hotter and less dynamically responsive at the time of the instability. That may be the reason why component C does not appear to be involved in the boxy/peanut structure of the bulge. Similarly, component D has a mean [Fe/H] of about −1.20, like that of the metal-weak thick disc near the Sun (e.g. Chiba & Beers 2000; Carollo et al. 2010).

Component E with mean [Fe/H] = −1.7 appears in small fractions at about the same [Fe/H] in most fields and shows no evidence for an abundance gradient with latitude. It is slowly rotating (see Ness et al., in preparation) and is probably an entirely separate population, not kinematically or chemically associated with the disc. We interpret the metal-poor component E stars as inner halo stars on orbits of high eccentricity passing through the bulge region. They may correspond to the high-eccentricity stars at an abundance of about −1.7 identified near the Sun by Chiba & Beers (2000). Simulations (e.g. Diemand, Madau & Moore 2005; Brook et al. 2007; Tumlinson 2010) show that the oldest stars in the Galaxy, formed in fragments from high-σ density fluctuations, should now inhabit the inner Galaxy. In this context, the stars of component E in the inner Galaxy would be worth further spectroscopic study at higher resolution.

We now return to the metal-rich component A, which is more prominent closer to the Galactic plane but at the same time shows most strongly the split red clump which we associate with the boxy/peanut structure of the bulge. The MDF of the solar neighbourhood may again provide a useful reference. In the solar neighbourhood, the youngest disc stars (ages < 2 Gyr) are relatively metal rich, have a relatively narrow range of abundance and are kinematically colder than the intermediate-age disc stars. The larger width of the MDF for the intermediate-age stars in the solar neighbourhood is currently attributed to the secular effects of radial migration, which radially smoothes the Galactic abundance gradient and so increases the width of the local MDF. In the same way, the more metal-rich stars of the early inner thin disc were probably kinematically colder and dynamically more responsive, thus preferentially mapped into boxy/bulge-supporting orbits.

Alternatively, it is possible that component A comes, via a second later and less vertically extended bar-forming instability event, from the more metal-rich thin disc which had continued to undergo star formation and chemical evolution since the first instability episode. Simulations have shown that a second buckling episode can occur (Athanassoula 2005b; Martinez-Valpuesta, Shlosman & Heller 2006), well after the first one, and that it leads to a vertically thicker structure than given by the first only buckling event. If we make the plausible assumption that component B forms earlier from the older metal-poor disc, and component A forms later from the metal-rich disc, then component B may well have had sufficient time to undergo a second buckling, while component A had not. This offers a possible explanation for why component B is more vertically extended and more metal poor than A. A secondary instability could also be triggered by a significant accretion event. Models with two discs of different scale heights may also offer more insight. The current N-body models we are working with do not show a secondary instability event and have a single disc only; we therefore cannot test these effects.

Component A attenuates in density by a factor of 5 between b = −5° and −10°. It is possible that the original bar formed in the early disc may have gradually entrained the chemically evolving thin disc during and after the instability process. In this way, the more metal-rich stars could have experienced only the later fraction of the bar-buckling process, and remained more closely bound to the Galactic plane.

In these alternative scenarios, we are associating component A with the inner thin disc which has continued to undergo star formation and chemical evolution. We can trace this component out to large Galactic radii, where it becomes gradually more metal poor (Table 5), reflecting the abundance gradient of the Galactic disc. In the central region however this component behaves in a clearly bar-like way, reflected in the split clump (see Section 6.2) and the kinematics as a function of longitude (Ness et al., in preparation).

7.2 The α-enhancement of the components

Fig. 18 shows that the [α/Fe]–[Fe/H] relation for stars in the inner Galaxy follows the pattern in the solar neighbourhood, except that the most metal-rich stars in the inner regions do show a weak α-enhancement. α-enhancement is a measure of the rate of chemical evolution. The more metal-poor boxy/peanut component B, which we associate with the older pre-bulge disc, is α-enhanced above the level of the more metal-rich bulge component A which we associate with the younger pre-bulge thin disc. This is consistent with our argument that the bulge stars of component B were formed earlier as part of the older disc. Their α-enrichment would be present not only in the bulge population but also for the other stars of the early disc which are still now in the disc.

The more metal-poor components in our survey, with [Fe/H] < −0.5, include the thick disc and metal-weak thick disc. Although there are several possibilities for the formation of the thick disc, the stars of components C and D may also have originated from the even earlier thin disc, heated kinematically by early minor merger activity before the instability event (e.g. Freeman 2008). They would be the oldest stars of the original disc, and their level of α-enhancement indicates that the duration of their star formation was short. Some of these more metal-poor stars may also have been puffed up into the bar/bulge by the instability event itself. The α-enhancement for all stars with [Fe/H] < −0.5 is fairly flat out to R_G = 30 kpc. As their α-enhancement is directly related to the duration of their star formation event, this suggests that the more metal-poor stars which contribute to the bulge/bar, the thick disc and metal-weak thick disc may have had a common origin in the oldest stars of the initial Galactic thin disc. In this picture, the different MDF components are fossil remnants from the dynamical redistribution of an initial disc system of the Milky Way. The abundances and α-enhancement are proxies for the time at which the stars of each component formed and the rate of their formation. Chemodynamical models may be able to link the abundances and time-scales for this formation and structural dispersion.

Fig. 18 shows the clearly bimodal distribution of abundances between the two boxy/peanut components A and B in our b = −5° fields. Component A has a mean [α/Fe] value of 0.16, and the old disc component B has a mean [α/Fe] of 0.23. Several authors have argued that the bulge is an α-enhanced metal-rich population (Fulbright, McWilliam & Rich 2007; Lecureur et al. 2007; Meléndez et al. 2008; Zoccali et al. 2008a). We see that the metal-rich component A is only weakly α-enhanced, and the thicker α-enhanced bulge component B is not so metal rich.

From Fig. 19, we can compare the α-enhancement of stars with similar [Fe/H] abundances in the bulge and in the disc at larger R_G. For stars with [Fe/H] > −1, which we would identify as originating from the disc, the mean [α/Fe] values are very similar in the bulge (Figs 19a and b) and the disc (Fig. 19c). For the more metal-poor stars with [Fe/H] < −1.0, the α-enhancement at large |R_G| (Fig. 19c) is similar to stars in the inner region |R_G| < 3.5 kpc. However, the dispersion in the α-enhancement is significantly lower in the inner bulge region (Fig. 19) compared to that at larger radii. This may reflect that the metal-poor population in the centre is more homogeneous than in the outer regions. For the very metal-poor stars with [Fe/H] < −2, the α-enhancement is significantly higher in the inner bulge region than in the disc at larger radii (see Fig. 17), indicating a more rapid star formation history for these very metal-poor stars in the inner regions of the Galaxy.

8 COMPARISON TO OTHER STUDIES

We have found an overall metallicity gradient in the whole MDF across our latitudes, comparable to that of Zoccali et al. (2008b). Our MDF at b = −10° and our α-enhancement results at our higher latitude also compare well with the results of Uttenthaler et al. (2012) for their high-resolution study at (l, b) = (0°, −10°). Note that Uttenthaler et al. (2012) also confirm the presence of the split clump only for the more metal-rich stars. Other spectroscopic studies of bulge stars have found evidence for multiple components in the MDF, and we now compare our results with theirs. Table 6 summarizes the mean abundances for the components found in the ARGOS survey compared with those from the work by Babusiaux et al. (2010), Bensby et al. (2010) and Hill et al. (2011). Each of these studies sees some but not all of the components found in the larger ARGOS sample.

8.1 Vertex deviations and the bulge

In the study of about 400 bulge K-giants and clump giants in Baade's window at b = −4° and a field near the minor axis at b = −6° by Babusiaux et al. (2010), our components A, B and D can be seen. They argue that both of the main scenarios of bulge formation have contributed to the bulge. Using proper motion and radial velocity data, they measure the vertex deviation of the velocity ellipsoid for their identified components. Their more metal-rich component, with [Fe/H] near 0.14 (our A), shows a significant vertex deviation and they conclude that this component is the bar/bulge that formed via the instability scenario. Their more metal-poor components (our B and D) show a smaller vertex deviation, and they associate these components with a spheroid and thick disc.

Babusiaux et al. (2010) also measured the radial velocity dispersion as a function of [Fe/H] for their fields at b = −4°, −6° and −12° near the minor axis of the bulge. They find that the dispersion increases with [Fe/H] in their b = −4° field (the more metal-rich stars are kinematically hotter) but decreases as a function of [Fe/H] in their higher latitude fields. Similarly, we find that at b = −5°, −7| ${.\!\!\!\!\!\!^{\circ}}$|5 and −10°, the velocity dispersion decreases with [Fe/H] along the minor axis (and indeed at all longitudes; see Ness et al., in preparation). This change with latitude of the sense of the abundance–(velocity dispersion) relationship indicates a structural change of the bulge populations at around b = −4°.

We now discuss the vertex deviation measurements of the bulge near the minor axis and compare these results to our N-body model. We do not have proper motions for our stars so cannot independently measure the vertex deviations. The vertex deviation of the stars in Baade's window has been measured by Soto, Rich & Kuijken (2007) and Babusiaux et al. (2010). Soto et al. (2007) conclude from their measurements that stars with [Fe/H] > −0.5 are supporting a bar. They find no significant vertex deviation for stars with [Fe/H] < −0.5. This makes sense in the context of our finding of the metallicity dependence of the split clump and assignment of components to populations. Soto et al. (2007) find a trend of decreasing velocity dispersion with increasing metallicity as do Babusiaux et al. (2010). Babusiaux et al. (2010) also find the largest vertex deviation for the more metal-rich stars which we would associate with component A. Their vertex deviation in the lr plane (denoted by |$l_{\text{v}}$|⁠) for stars with [Fe/H] < −0.14 is −13° ± 9°. For the more metal-rich stars with [Fe/H] > −0.1, the vertex deviation is about −40°. They argue that the most metal-rich stars represent the bar/bulge in the central region, and the stars with lower metallicities are associated with an old spheroid.

We consider first the vertex deviation of stars belonging to our boxy/peanut bulge component B, with 0 > [Fe/H] > −0.5. The |$l_{\text{v}}$| vertex deviation measured by Babusiaux et al. (2010) for stars in this abundance range is about −30°. This agrees well with the estimate from our bar/bulge model, for which we find a vertex deviation |$l_{\text{v}} = -26^\circ \pm 2^\circ$| for all stars within a galactocentric radius of 0.5 kpc.

Now we address the vertex deviation of stars belonging to our boxy/peanut bulge component A, with [Fe/H] between 0 and 0.5. For stars in this abundance range, Babusiaux et al. (2010) measure a vertex deviation in Baade's window of |$l_{\text{v}} = -43^\circ \pm 5^\circ$|⁠. We examined six models with different shapes and scale heights of the central bulge, scaled to match the kinematics of the Milky Way, and found the vertex deviation near the minor axis to be fairly insensitive to the scale heights of the bulge/bar and also to the height above the plane. To reproduce this higher |$l_{\text{v}}$| in our models, we needed to rotate the bulge/bar in the model from the previously adopted 20° to 45° with respect to the line of sight. The vertex deviation in the model for a bar angle of 45° is |$l_{\text{v}} = -45^\circ \pm 2^\circ$| similar to that derived by Babusiaux et al. (2010). This is an interesting result in light of the work of Cabrera-Lavers et al. (2008) and others (e.g. Weinberg 1992; Hammersley et al. 1994) who argue for the existence of a long bar in the Milky Way from analysis of 2MASS data (Skrutskie et al. 2006). Using red clump giants as tracers, Cabrera-Lavers et al. (2008) present evidence for a long bar with a semilength of 4 kpc, ending around l = 28°. This long bar inclined at 45° coexists with the thicker bulge/bar inclined at 20° and ≈3.5 kpc in length (Gerhard 2002) (our component B). It has been previously suggested (e.g. Cabrera-Lavers et al. 2008) that the inner Galaxy is dominated at higher latitudes by the bulge and at lower latitudes by the long bar.

We are not suggesting that the vertex deviation in our model is strong evidence for the coexistence of a boxy/peanut bulge made of stars with 0 > [Fe/H] > −0.5 and a long flatter bar that contains the metal-rich stars with [Fe/H] > 0 in the bulge. It is however an interpretation in light of the dependence of vertex deviation of the model on the bar orientation. Changing the angle however may be one of several ways to affect the vertex deviation. Given the magnitude distribution of the split clump at the higher latitudes as a function of [Fe/H] and the changing relationship between dispersion and [Fe/H] at b = −5°, it seems clear that there are two components to the bulge population. It is possible that these results could be explained by a secondary structure which is dominant at low latitudes as per Cabrera-Lavers et al. (2008). The metal-rich stars may be associated with a thin disc population associated with a bar orientated at 45° to the Sun–centre line.

8.2 Microlensed dwarfs in the bulge

A potential problem with our interpretation of component B as the primary component of the bulge comes from Bensby et al. (2011) who compiled abundances for 26 microlensed dwarfs in the direction of the bulge. The [Fe/H] values for their stars extend from −0.72 to +0.54. The MDF of their microlensed dwarfs shows two well-separated components corresponding to our components A and C (which we argue are the thin boxy/peanut bulge and the old thick disc). As for our sample, they find that the stars of component C are α-enhanced, while the metal-rich stars of component A have [α/Fe] ≈ 0.1, close to solar. They estimate isochrone ages for these stars and find that the component C stars are almost all old, around 10–12 Gyr, while the component A stars cover a wide age range from 2 to 14 Gyr. Their ages are consistent with the identifications of component A with the thin disc in a boxy/peanut bar structure, which has continued to undergo star formation, and with C being the old thick disc. Bensby et al., however, find only a small fraction of stars in the [Fe/H] range of what we believe to be the main bar/bulge (component B), around [Fe/H] = −0.25.

We note that the Bensby et al. (2011) stars are mostly closer to the plane than our lowest fields at b = −5°. At b = −5°, at a radial distance of ±0.5 kpc from the centre of the bulge, about 13 per cent of our stars have [Fe/H] between −0.2 and 0, but there is only one microlensed dwarf out of 26 in this abundance range. Along the minor axis at b = −5°, 40 per cent of our stars have [Fe/H] > 0, somewhat lower than the 50 per cent of microlensed stars in the same abundance interval which lie closer to the plane. This increase in the fraction of metal-rich stars closer to the plane is consistent with our finding that the fraction of the most metal-rich component A increases towards the plane. The number of stars in the metallicity range of component B (0 > [Fe/H] > −0.5) increases by a factor of only 1.2 from b = −10° to −5°, whereas the number of stars in the metallicity range of component A ([Fe/H] > 0) increases by a factor of 5.6. Furthermore, if A is associated with the thin disc, then there may be an increasingly mixed-age population of metal-rich stars at low latitudes that is not seen at higher latitudes (b < −5°).

The missing population B at the latitudes measured by Bensby et al. (2011) suggests that the component B stars close to the plane are greatly outnumbered by the metal-rich stars of component A. The region near the plane may be dominated by the flatter bar and host a younger, thin disc population which may be part of the bar. The low fraction of stars from component B in the Bensby et al. (2011) data suggests that the stars of the early thin disc with [Fe/H] ≈ −0.25 were mostly displaced to z-heights > 0.7 kpc during the boxy/peanut bulge buckling process. As a check on these speculations, it would be interesting to know the abundance distribution of microlensed dwarfs at higher latitudes, closer to ours. The probability of finding microlensed stars decreases rapidly, however, with latitude. It is also possible that the Bensby et al. (2011) stars do not lie in the bulge itself: the geometry of microlensing would favour their preferential location on the far side of the peanut/boxy bulge. In this case, one would expect to find a relatively few stars of component B (the dominant peanut/boxy bulge population) in the microlensed sample. To check this suggestion, we would need a realistic calculation of the dependence of microlensing optical depth on the position of stars in the inner regions of the Galaxy.

9 CONCLUSIONS

By measuring stellar parameters for the 25 500 stars in our survey, we are able to identify the nearby dwarfs and subgiants, determine distances to the giants and make a distance cut to separate out approximately 14 150 bulge candidates with |R_G| ≤ 3.5 kpc. The MDF of these 14 150 stars shows four to five distinct components (denoted by A–E). We find that the relative contributions of these components change with position in the bulge. Our data extend into the disc beyond the bulge, which helps to identify the nature of the MDF components and interpret how they originated from the early disc of the Milky Way.

The components show that the stars of the inner Galaxy are a composite population. The inner region of the Milky Way comprises stars with [Fe/H] > −0.5, which appear to belong to the boxy/peanut bulge, and stars with [Fe/H] < −0.5 which lie in the inner regions but are not part of the boxy/peanut structure. There are two spatially and chemically separated bulge components with [Fe/H] > −0.5 – (i) component A: a relatively thin boxy/peanut bulge with mean [Fe/H] ≈ 0.15 which we associate with the thin disc and (ii) component B: a thicker boxy/peanut bulge population with mean [Fe/H] ≈ −0.25.

The duality of the bulge is clear but the mechanism by which this duality has arisen is not so clear. N-body models of disc evolution show that the thick boxy/peanut component B is likely to have formed in a major instability event which occurred many Gyr ago. The thinner boxy/peanut bulge component A appears also to be associated with such an instability event. It may have originated in the same event from the initial stratification of the thin disc, or it might have formed after the thicker structure.

Chemically we associate the thin boxy/peanut bulge (component A) with the younger thin disc, and this component can be followed out to large Galactic radii beyond the bulge at low |z|-heights from the plane where it merges with component B. It shows a radial metallicity gradient, becoming more metal poor in the outer regions. Its gradient on the near side of the bulge from R_G =0 to 6 kpc is −0.02 ± 0.01 dex kpc⁻¹. Spatially it seems likely that component A comes from the redistribution of the younger thin disc in the inner regions into a boxy/peanut structure, via the bar-forming and bar-buckling instabilities of the disc. Our MDFs indicate that nearer the Galactic plane this population will be very dominant. It may also be mixed in age, due to ongoing star formation, with younger stars located nearest to the plane.

Component B is the primary bulge component that is more vertically extended than the metal-rich component A. Outside the bulge at high latitudes, where a very few metal-rich stars are present even inside the bulge, it represents the most metal-rich fraction of stars observed in our sample at large |z|-heights from the plane. It is associated with the early disc and from the instability scenario represents the primary population which formed the boxy/peanut bulge.

The distinct mean [Fe/H] and [α/Fe] values for components A and B suggest that they formed in two distinct star-forming episodes. Comparison of vertex deviations observed in the inner Galaxy and derived for the N-body models indicates that the thin boxy/peanut bar component A and the thick boxy/peanut bulge component B may lie at different angles to the Sun–centre line: 45° and 20°, respectively. The two components A and B show the same vertical abundance gradient of −0.08 ± 0.05 and −0.08 ± 0.04 dex kpc⁻¹, respectively. An even smaller vertical abundance gradient of −0.04 ± 0.02 dex kpc⁻¹ is seen for the thick disc component C in the inner Galaxy. These small vertical gradients need not be interpreted as favouring a merger scenario for the formation of the bulge, or formation via dissipational collapse. They are likely to be residual gradients of the disc surviving from the instability event. The survival of vertical gradients in the initial disc needs to be confirmed by a more detailed study of the way in which the instability process in the N-body models maps the initial phase space into the post-instability phase space. We argue that the much larger overall vertical abundance gradient in the bulge region, of about −0.45 to −0.6 dex kpc⁻¹, comes primarily from changing contributions of components A–C with height above the plane (see Fig. 11).

We find that the bulge is not a single population. Its metallicity components reflect the Galactic components which were present before the boxy/peanut bulge formation, and the properties of these components reflect the dynamical processes that altered the spatial distribution of stars primarily in the early disc. Stars with [Fe/H] > −0.5 are involved in the boxy/peanut bulge structure, and stars with [Fe/H] < −0.5 are present in the inner bulge regions and contribute to the total population of the inner Galaxy, particularly at higher latitudes. Our thin and thick boxy/peanut bulge/bar structures are both similar to the boxy/peanut bar/bulge seen in our instability-driven N-body model, but with different scale heights and correspondingly different kinematics. Our more metal-poor components C, D and E, with [Fe/H] < −0.5, can be identified with the pre-instability thick disc, metal-weak thick disc and halo populations.

The presence of the two kinematically and chemically distinct components A and B in the boxy/peanut bulge is evidence that internal evolution is the dominant process in the formation of the bulge. It seems unlikely that major mergers have been an important part of the formation of the bulge of the Milky Way. Such disruption to the inner region would probably smooth the two bulge components into a continuous population. The survival of the two components of the bulge constrains the level of any major disturbance that could have occurred in the inner Galaxy after the disc instability event. At this stage, we do not attempt to place limits on the contribution of any merger-generated ‘classical’ bulge component. Further N-body modelling is needed to quantify this issue.

We thank the staff at the Australian Astronomical Observatory, who have made this project possible. This publication makes use of data products from the Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation. This work has been supported by the RSAA and Australian Research Council grant DP0988751. EA gratefully acknowledges financial support by the CNES and by the European Commission through the DAGAL Network (PITN-GA-2011-289313). JBH is supported by an ARC Federation Fellowship. GFL thanks the Australian Research Council for support through his Future Fellowship (FT100100268) and Discovery Project (DP110100678). RRL gratefully acknowledges financial support from FONDECYT, project No. 3130403. LLK is supported by the Lendület programme of the Hungarian Academy of Sciences and the Hungarian OTKA Grants K76816, MB08C 81013 and K83790.

REFERENCES

APPENDIX A:

Regarding our red clump magnitude calibration: the Alves (2000) calibration is for stars in the range 0.0 > [Fe/H] > 0.5, which is the metallicity range of the majority of our stars in the bulge. Alves investigated the metallicity dependence on the absolute magnitude of the clump and recommended no [Fe/H] correction be made for the stars in this abundance range. Laney, Joner & Pietrzyński (2012) similarly reported no dependence on metallicity and an absolute red clump magnitude of M_K = −1.61, in agreement with Alves (2000).

We have adopted his M_K over a broader metallicity range, because there is disagreement in the literature about the existence of a K-band magnitude–metallicity dependence. It is consistently reported to be weak at most, and there is disagreement about the factor in the equation for the metallicity dependence. We believe that it is appropriate to adopt a mean value of −1.61 for calculating clump distances for our stars (see Valentini & Munari 2010 for a review). Even if a metallicity dependence exists, its effect on our distances is likely to be small. For example, with Udalski's I-band calibration of a 0.13 mag dex⁻¹ dependence on abundance, the distance of a star with [Fe/H] = −1.0 would change by about 6 per cent. We note that 95 per cent of our stars have [Fe/H] > −1.0.

Pietrzyński, Gieren & Udalski (2003) looked at the red clump magnitude dependence on age and [Fe/H] in J and K and found only a low or no dependence of the K magnitude on [Fe/H] over the wide metallicity range of the Large Magellanic Cloud, Small Magellanic Cloud and the Carina and Fornax dwarf galaxies. They concluded that the mean K-band magnitude of the red clump is an appropriate distance indicator for a wide range in [Fe/H] and age. van Helshoecht & Groenewegen (2007) similarly found that the absolute magnitude is an appropriate distance indicator in the range 0.4 > [Fe/H] > 0.5, although they adopted a slightly fainter value of M_K = −1.57 for the mean absolute magnitude of the clump.

We note that theory predicts a metallicity dependence of the absolute magnitude of the clump (i.e. Girardi & Salaris 2001) and this has been supported by earlier observations (i.e. Sarajedini 1999; Twarog, Anthony-Twarog & Bricker 1999). Udalski (1998) found a weak metallicity dependence in the I band in the range 0.6 < [Fe/H] < 0.2 using red clump stars in 15 Magellanic cloud clusters. However, he concluded that the absolute I-band magnitude of the red clump is independent of age for ages 2–10 Gyr. In summary, we will adopt M_K = −1.61 for the clump, independent of age and abundance.