Testing the Universality of the Stellar IMF with Chandra and HST

Understanding the stellar initial mass function (IMF) has major implications for a variety of astrophysical problems. The IMF plays a central role in converting the observed properties of galaxies (e.g., luminosity and color) into physically meaningful quantities like stellar mass and star formation rate. An assumed universality of the IMF has contributed to our current observational-based paradigm of how galaxies formed and evolved throughout the history of the universe, what fraction of the universe's mass is tied up in stellar baryons, and the number of compact objects in the universe. A clear understanding of the form and universality of the IMF in external galaxies is therefore a central goal in modern astrophysics (see, e.g., Bastian et al. 2010 for a review).

In recent years, there has been mounting evidence indicating that the stellar IMFs for elliptical galaxies vary with galaxy mass. For example, Cappellari et al. (2012) used detailed dynamical models and two-dimensional stellar kinematic maps of 260 early-type galaxies from the ATLAS${}^{{\rm{3D}}}$ project to show that the mass-to-light ratios (M/L) for elliptical galaxies increase with increasing stellar velocity dispersion, σ, consistent with a scenario where the IMF changes with galaxy mass. Such a finding has been independently noted in galaxy lensing measurements of M/L for a variety of σ (e.g., Auger et al. 2010; Treu 2010). Specifically, these findings indicate that relatively low-mass, early-type galaxies (σ ≲ 100 km s⁻¹) have M/L values consistent with standard Milky Way-like IMFs (Kroupa 2001; Chabrier 2003). However, relatively massive ellipticals (σ ≈ 300 km s⁻¹) have M/L that are larger than those predicted from standard IMFs, and can be consistent with either "bottom-heavy" IMFs ($\alpha \approx 2.8$, yielding more low-mass stars with higher M/L) or "top-heavy" IMFs ($\alpha \approx 1.5$, yielding more remnants with lower luminosities; Cappellari et al. 2012).

van Dokkum & Conroy (2010, 2011, 2012) and Conroy & van Dokkum (2012) showed that the spectra of massive ellipticals have strong Na i and Ca ii absorption features that are indicative of a large population of very low-mass stars ($\lesssim$0.3 ${M}_{\odot }$; Saglia et al. 2002), favoring the bottom-heavy IMF interpretation for these galaxies. These features had been known for quite some time (see, e.g., Hardy & Couture 1988; Cenarro et al. 2003), but higher-quality data and spectral stellar synthesis modeling in recent years have given more confidence to the interpretation that the stellar IMFs in these galaxies are likely to be bottom-heavy and consistent with $\alpha \approx 2.8$.

Despite the above evidence for a bottom-heavy IMF in massive ellipticals, some inconsistencies remain. For example, Smith & Lucey (2013) found that a bottom-heavy IMF is incompatible with the stellar mass obtained via strong lensing for a galaxy with $\sigma \approx 300$ km s⁻¹. Smith (2014) has also shown, intriguingly, that there is no correlation between results based on near-IR indices and those based on dynamics, on a galaxy-by-galaxy basis. Additionally, Weidner et al. (2013) showed that a time-invariant and bottom-heavy IMF was incompatible with the observed chemical enrichment of massive ellipticals, and underpredicted the number of stellar remnants (in the form of LMXBs) observed in globular clusters (GCs). To resolve this, they proposed a time-dependent IMF model that would evolve from a top-heavy to bottom-heavy form. In this scenario, the star formation histories of massive ellipticals would have included early phases of intense starburst activity that produced many massive stars (and a top-heavy IMF), which chemically enriched the galaxies and created turbulence in their interstellar mediums. Subsequent star formation within these environments would be more greatly fragmented and lead to the formation of preferentially smaller stars and bottom-heavy IMFs.

Other authors have taken somewhat different approaches to resolving the above incompatibility issue. Ferreras et al. (2015) focused on the functional form of the IMF and tested a model where the low-mass and high-mass slopes of a broken power-law IMF were independently varied. They found that regardless of the IMF slope parameters, a time-independent IMF (independent of galaxy mass) could not simultaneously reproduce the observed chemical enrichment and gravity-sensitive spectral features in massive ellipticals. Finally, Martín-Navarro et al. (2015) explored a radial and time-dependent form of the IMF. They argued that in massive ellipticals ($\sigma \approx 300$ km s⁻¹), the radial trends in chemical enrichment and gravity-sensitive spectral features could be well described by a time-dependent IMF, like that proposed by Weidner et al. (2013), in the central regions of these galaxies. In contrast, a Milky Way-like IMF was proposed in the galactic outskirts, based on the lack of these spectral signatures relative to the galactic interior. The IMF was then claimed to be a local phenomenon, reflecting a two-stage formation history for these galaxies. Other evidence for radial trends in spectral features follows similar rationales (La Barbera et al. 2016a, 2016b).

The studies above have helped to constrain variations at the low-mass end of the elliptical-galaxy IMF ($\lesssim$0.5 ${M}_{\odot }$) with galaxy mass but do not place strong constraints on variations at the high-mass end of the IMF. To trace the high-mass end of the IMF requires studying the most massive stars in a galaxy; however, these stars evolve quickly into black holes (BHs) and neutron stars (NSs), making them difficult to observe. Fortunately these stellar remnants are found in binary star systems, where a lower-mass and longer-lived companion star eventually becomes a "donor" star, losing mass to the compact object. The system itself radiates nearly all of its energy in the form of X-rays, making these previously undetectable high-mass remnants available to direct observation. Therefore, an efficient and effective way to constrain variations at the high-mass end of the IMF ($\gtrsim$8 ${M}_{\odot }$) with galaxy mass is to test for corresponding variations in the prevalence of these low-mass X-ray binary (LMXB) systems in elliptical galaxies. With this in mind, Peacock et al. (2014, hereafter, P14) calculated (using population synthesis calculations from Maraston 2005) that a bottom-heavy single power-law IMF, with $\alpha \approx 2.8$, would yield a factor of ≈3 times fewer LMXBs per unit K-band luminosity ($N/{L}_{K}$) than a standard Kroupa (2001) IMF, implying a dramatic decline in $N/{L}_{K}$ with increasing σ. Using archival Chandra X-ray Observatory (Chandra) and Hubble Space Telescope (HST) data on eight nearby elliptical galaxies, P14 found instead that $N/{L}_{K}$ was constant with their luminosity and velocity dispersion; however, their test included only one low-mass elliptical galaxy at σ < 150 km s⁻¹, where $N/{L}_{K}$ was expected to be largest.

In this paper, we augment the P14 study of eight elliptical galaxies by adding Chandra constraints for five new low-mass ellipticals (σ = 78–110  km s⁻¹) that are predicted to have standard IMFs that differ from those of the massive ellipticals (see, e.g., Cappellari et al. 2012). With this expanded sample of 13 galaxies, we are able to place robust constraints on how $N/{L}_{K}$ varies across the full range of velocity dispersion (σ = 80–300 km s⁻¹) where significant variations in $N/{L}_{K}$ are plausibly expected.

Our sample of low-mass elliptical galaxies was derived from the samples of Halliday et al. (2001) and Cappellari et al. (2006), which contain velocity dispersion data for several nearby ellipticals. We limited our sample to galaxies with D < 20 Mpc, so that we could easily resolve LMXB populations, and σ ≲ 110 km s⁻¹, to focus on the galaxies that are likely to have standard IMFs that differ from those of the already well-studied high-mass ellipticals. In order to guard against potential contamination from young high-mass X-ray binaries (HMXBs), we further restricted our sample to galaxies that had negligible signatures of star formation rate (SFR < 0.001 ${M}_{\odot }$ yr⁻¹), as measured by 24 μm Spitzer data (Temi et al. 2009). Our final sample of six low-mass elliptical galaxies is summarized in Table 1.

Table 1.  Properties of the Sample of Low-mass Elliptical Galaxies

	D	σ	a	b		${N}_{{\rm{H}}}$	HST ACS Data		${t}_{\exp }$	${N}_{\mathrm{LMXB}}$	${N}_{{\rm{X}},\mathrm{GC}}$	${N}_{{\rm{X}},\mathrm{bkg}}$
Source Name	(Mpc)	(km s−1)	(arcmin)		log(${L}_{K}/{L}_{K,\odot }$)	(1020 cm−2)	(Blue Filter)	(Red Filter)	(ks)	(Field)	(GCs)	(Background)
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)	(10)	(11)	(12)	(13)
NGC 4339	16.0	100.0	1.3	1.1	10.3a	1.62	F606W	⋯	33.6	2a	0a	1a
NGC 4387	17.9	97.0	0.9	0.6	10.2	2.73	F475W	F850LP	38.7	1	0	1
NGC 4458	16.4	85.0	0.9	0.7	10.0	2.63	F475W	F850LP	34.5	2	0	1
NGC 4550	15.5	110.0	1.3	0.5	10.2	2.60	F475W	F850LP	25.8	6	1	1
NGC 4551	16.1	95.0	1.1	0.7	10.2	2.59	F475W	F850LP	26.6	0	1	0
NGC 7457b	12.9	78.0	2.6	1.4	10.3	5.49	F475W	F850LP	37.7	1	0	2
Total Sample										12	2	6

Notes. Column (1): Target galaxy name. Column (2): Distance as given by Cappellari et al. (2013). Column (3): Velocity dispersion from either Halliday et al. (2001) or Cappellari et al. (2006). Columns (4) and (5): 2MASS-based K-band major and minor axes of the galaxy from Jarrett et al. (2003). Column (6): Logarithm of the K-band luminosity. Column (7): Galactic column density of neutral hydrogen. Columns (8) and (9): The available HST imaging (via either ACS or WFPC2 imaging) for "blue" and "red" filters, respectively, which were used to identify and characterize optical counterparts to X-ray sources. Column (10): Total Chandra exposure time. All galaxies were imaged using ACIS-S and had coverage over the entire galactic extents as defined in Columns (4) and (5). Columns (11)–(13): The number of field LMXBs (Column (11)), GC LMXBs (Column (12)), and background or central AGN candidates (Column (13)) within the extent of each galaxy (as defined in Columns (4) and (5)) that had 0.5–7 keV fluxes exceeding that of a source with ${L}_{{\rm{X}}}={10}^{38}$ erg s⁻¹ at the distance provided in Column (2). X-ray sources were classified using HST data and the procedure outlined in Section 2.2 of P14.

^aFor NGC 4339, we considered 74% of the reported K-band luminosity that corresponded to the region of the galaxy covered by the HST WFPC2 footprint detailed in Section 2. For Columns (11)–(13), we considered only X-ray sources detected in this covered region. ^bObservations of NGC 7457 were represented in both the sample of P14 and this study. However, due to the availability of new HST data for NGC 7457 in this study, results based on our analysis were used throughout this paper.

Download table as:  ASCIITypeset image

We conducted new Chandra observations for four of the six galaxies to reach 0.5–7 keV point-source detection limits of ${L}_{{\rm{X}}}\approx {10}^{38}$ erg s⁻¹ (chosen for consistency with the sample from P14), after combining with archival Chandra data. All observations were conducted using ACIS-S, which covers the entire K-band-defined areal footprints of all of the galaxies in our sample (see Table 1 and Jarrett et al. 2003). The only exception was NGC 4339, which had only partial HST coverage. For this galaxy, we utilized only the background-subtracted K-band flux within the overlapping HST footprint, and only LMXBs detected in this region were included in our analyses. This resulted in our using only 74% of the K-band luminosity reported in Column (6) when computing $N/{L}_{K}$ for this galaxy.

Data reduction for our sample of galaxies closely follows the procedure outlined in Sections 2.1 and 2.2 of Lehmer et al. (2013), with our reduction being performed with CIAO v. 4.7 and CALDB v. 4.6.7. We reprocessed events lists from level 1 to level 2 using the script chandra_repro, which identifies and removes events from bad pixels and columns, and filters events lists to include only good time intervals without significant flares and non-cosmic-ray events corresponding to the standard ASCA grade set (grades 0, 2, 3, 4, 6). For galaxies with more than one observation, we combined events lists using the script merge_obs. We constructed images in three X-ray bands: 0.5–2 keV, 2–7 keV, and 0.5–7 keV. Using our 0.5–7 keV images, we utilized wavdetect at a false-positive probability threshold of 10⁻⁵ to create point-source catalogs. We converted 0.5–7 keV point-source count rates to fluxes assuming an absorbed power-law spectrum with a photon index of ${\rm{\Gamma }}=1.5$ and Galactic extinction (see Column (7) in Table 1). We treated the hot gas component as negligible within the small area of the sources, because these galaxies have very little diffuse emission in total, typical of low-mass ellipticals like those studied here (e.g., O'Sullivan et al. 2001). Our choice of photon index reproduces well the mean 2–7 keV to 0.5–2 keV count-rate ratio of the detected point sources in our sample, a value calculated using stacking analyses (see, e.g., Lehmer et al. 2016 for details).

The total number of sources within each galaxy footprint that had fluxes brighter than the limit ${L}_{{\rm{X}}}={10}^{38}$ erg s⁻¹ spans the range 1–8 (see Columns (11)–(13) of Table 1). Our test of the variation in the IMF with galaxy mass is sensitive to the field LMXB population that forms within the galactic stellar populations. To search for potential contaminants from unrelated X-ray-detected objects, including, e.g., background active galactic nuclei (AGNs), foreground Galactic stars, central low-luminosity AGNs, and LMXBs formed through dynamical processes in GCs, we made use of archival HST data (Columns (8) and (9) in Table 1). Contaminating source populations were identified following the procedure outlined in Section 2.2 of P14, which makes use of HST source colors and morphologies to identify and classify counterparts. Regardless of their nature, we rejected X-ray sources with optical counterparts from further consideration. In total, we rejected eight of the X-ray sources within the galactic footprints (see Columns (12) and (13) of Table 1) that had optical counterparts, which left 12 candidate field LMXBs with ${L}_{{\rm{X}}}\gt {10}^{38}$ erg s⁻¹ in our sample.

In Figure 1(a), we show the number of field LMXBs with ${L}_{{\rm{X}}}\gt {10}^{38}$ erg s⁻¹ per unit K-band luminosity, ${N}_{38}/{L}_{K}$, versus velocity dispersion, σ, for our sample of low-mass ellipticals combined with the P14 sample. The "variable" and "invariant" IMF models from P14 are shown in Figure 1(a) with blue and red curves, respectively. The variable model assumes that the IMF makes a transition from Kroupa at $\sigma =95$ km s⁻¹ (i.e., $\mathrm{log}\sigma =1.98$) to a single power law with slope $\alpha =2.8$ at $\sigma =300$ km s⁻¹ (i.e., $\mathrm{log}\sigma =2.48;$ see details below), while the invariant model assumes a Kroupa IMF at all σ.

Zoom In Zoom Out Reset image size

At first inspection, it seems that the specific frequency of LMXBs indicates that the IMFs of elliptical galaxies are consistent with a single "universal" IMF. However, it is important to note that the LMXB population is tracing only the remnant population from stars with m ≳ 8 ${M}_{\odot }$. All that can be reliably inferred from the apparent constancy of ${N}_{38}/{L}_{K}$ across σ is that the single power-law slope of the IMF for stars of $m\gt 0.5$ ${M}_{\odot }$ does not undergo strong variations for ellipticals of all velocity dispersions. In fact, it may be possible that the low-mass end does vary strongly with velocity dispersion, as found in the literature (see references in Section 1).

With this in mind, we sought to extend P14 by exploring the space of acceptable parameters that could constrain how the IMF might vary, by constructing a suite of IMF models and comparing their predicted ${N}_{38}/{L}_{K}$ versus σ tracks with our constraints shown in Figure 1(a). In this process, we followed a similar, but more generalized procedure to that outlined in Section 4.2 of P14. Below, we summarize our model.

We began by considering a broken power-law form for an IMF corresponding to high-mass ellipticals:

$$\begin{eqnarray}\displaystyle \frac{{dN}}{{dm}}={N}_{0}\left\{\begin{array}{ll}{2}^{({\alpha }_{2}-{\alpha }_{1})} & {m}^{-{\alpha }_{1}}\,\quad 0.1\,{M}_{\odot }\lt m\lt 0.5\,{M}_{\odot }\\ & {m}^{-{\alpha }_{2}}\,\quad m\gt 0.5\,{M}_{\odot },\end{array}\right.\end{eqnarray} \tag{ 1 }$$

where ${\alpha }_{1}=1.3$ and ${\alpha }_{2}=2.3$ for a Kroupa IMF. N₀ is a constant of normalization, which, in our procedure, normalizes the 0.1–100 ${M}_{\odot }$ integrated IMF to 1 ${M}_{\odot }$. Therefore, N₀ varies with ${\alpha }_{1}$ and ${\alpha }_{2}$.

Next, we constructed a grid of IMFs over ${\alpha }_{1}=1$–5 and ${\alpha }_{2}=1$–3.5 with 801 and 501 steps of 0.005, respectively, thus resulting in a grid of n = 401,301 unique IMFs. For the ith IMF, we quantified the K-band mass-to-light ratio, ${(M/{L}_{K})}_{i}$, by running the stellar population synthesis code PÉGASE (Fioc & Rocca-Volmerange 1997, 1999), adopting for consistency the assumptions in P14 of a star formation history comprising a single burst of age 10 Gyr and solar metallicity. In this procedure, M represents the total initial stellar mass and therefore does not vary with age. We defined the ratio of the ith mass-to-light ratio to that of the Kroupa case as

$$\begin{eqnarray}&&{R}_{({\text{}}M/L),i}=\displaystyle \frac{{(M/{L}_{K})}_{i}}{{(M/{L}_{K})}_{{\rm{kro}}}}=\displaystyle \frac{{L}_{K,\mathrm{kro}}}{{L}_{K,i}}.\end{eqnarray} \tag{ 2 }$$

In order to keep ${R}_{({\text{}}M/L)}$ consistent with observations of how the dynamical and stellar M/L varies with σ, we required that ${R}_{({\text{}}M/L),i}\lt 3.0$ (e.g., Cappellari et al. 2012; Conroy & van Dokkum 2012). This requirement on ${R}_{({\text{}}M/L)}$ resulted in the discarding of models with specific combinations of ${\alpha }_{1}$ and ${\alpha }_{2}$ (notably, models with ${\alpha }_{1}\geqslant 3.9$ or ${\alpha }_{2}\leqslant 1.6$ were rejected; see Figure 1(d)). This constraint limited the number of models considered to a total of n = 198,298, which we utilize hereafter.

Because the prevalence of the LMXB population is sensitive to the underlying population of compact objects, i.e., NSs and BHs, which are remnants of >8 ${M}_{\odot }$ stars, we calculated the number of stars per solar mass that become compact objects for a given IMF by integrating the IMF from 8 to 100 ${M}_{\odot }$:

$$\begin{eqnarray}&&{N}_{\mathrm{CO},i}={N}_{0,i}{\int }_{8}^{100}{m}^{-{\alpha }_{2,i}}{dm}.\end{eqnarray} \tag{ 3 }$$

where ${N}_{0,i}$ is the ith normalization factor. These mass fractions allow us to compute the K-band luminosity-normalized ratio of expected NSs and BHs generated by the ith IMF to that of a Kroupa IMF from

$$\begin{eqnarray}&&{R}_{{\rm{CO}},i}\equiv \displaystyle \frac{{({N}_{\mathrm{CO}}/{L}_{K})}_{i}}{{({N}_{{\rm{CO}}}/{L}_{K})}_{\mathrm{kro}}}=\displaystyle \frac{{N}_{\mathrm{CO},i}}{{N}_{{\rm{CO,}}{\rm{kro}}}}{R}_{({\text{}}M/L),i}.\end{eqnarray} \tag{ 4 }$$

From these quantities, we construct a variable IMF model that varies smoothly with σ, creating a bridge from the Kroupa IMF of the low-mass ellipticals to the IMF of the high-mass ellipticals (i.e., the ith IMF). We define our variable IMF function over the range σ = 95–300 km s⁻¹, within which we require that the IMF varies as a function of σ following the relation ${({\text{}}M/L)}_{\sigma }\propto {\sigma }^{0.72}$ of Cappellari et al. (2013). Under these assumptions, we can quantify the fraction of the variable IMF that is composed of the ith IMF as a function of σ to reflect the desired result that galaxies with $\sigma \ \leqslant$ 95 km s⁻¹ have no contribution from the ith IMF, but will be composed of only a Kroupa IMF, while galaxies with $\sigma \ \geqslant$ 300 km s⁻¹ will be composed of solely the ith IMF component:

$$\begin{eqnarray}F(\sigma )=\left\{\begin{array}{ll}0 & (\sigma \leqslant 95\,\mathrm{km}\,{{\rm{s}}}^{-1})\\ \tfrac{{\sigma }^{0.72}-{95}^{0.72}}{{300}^{0.72}-{95}^{0.72}} & (95\,\mathrm{km}\,{{\rm{s}}}^{-1}\lt \sigma \lt 300\,\mathrm{km}\,{{\rm{s}}}^{-1})\\ 1 & (\sigma \geqslant 300\,\mathrm{km}\,{{\rm{s}}}^{-1}),\end{array}\right.\end{eqnarray} \tag{ 5 }$$

We define the complementary fraction ${F}_{\mathrm{kro}}$,

$$\begin{eqnarray}&&{F}_{{\rm{kro}}}(\sigma )=1-F(\sigma ),\end{eqnarray} \tag{ 6 }$$

which is then the fraction of the variable IMF that is composed of a Kroupa IMF.

We combine Equations (4)–(6) and define a composite function:

$$\begin{eqnarray}{R}_{{\rm{comp}},i}(\sigma ) & \equiv & \displaystyle \frac{{({N}_{\mathrm{CO}}/{L}_{K})}_{\sigma ,i}}{{({N}_{{\rm{CO}}}/{L}_{K})}_{\mathrm{kro}}}={R}_{{\rm{CO}},i}F(\sigma )+{F}_{{\rm{kro}}}(\sigma )\\ & = & 1-(1-{R}_{{\rm{CO}},i})F(\sigma ),\end{eqnarray} \tag{ 7 }$$

which represents the σ-dependent number of compact objects per unit K-band luminosity compared to the Kroupa IMF. We then arrive at a function that can be used to predict the number of observed LMXBs per unit K-band light:

$$\begin{eqnarray}{\left(\displaystyle \frac{{N}_{{\rm{X}}}}{{L}_{K}}\right)}_{i} & = & {\xi }_{i}\left(\displaystyle \frac{{N}_{{\rm{CO,}}{\rm{kro}}}}{{L}_{K,{\rm{kro}}}}\right){R}_{{\rm{comp}},i}(\sigma )\\ & = & {\gamma }_{{\rm{X}},i}{R}_{{\rm{comp}},i}(\sigma ).\end{eqnarray} \tag{ 8 }$$

Here ${\xi }_{i}$ represents the luminosity-dependent fraction of the population of compact objects that is actively involved in an LMXB phase, derived using the ith model. ${\gamma }_{{\rm{X}},i}$ is a fitted scaling factor for the ith model that allows us to express this quantity as an observed frequency of LMXBs by correcting for a portion of N_CO that are not LMXBs (e.g., BH–BH pairs, BH–NS pairs, etc.). We note that we assume that ξ is independent of σ, which simplifies the computation of our IMF models, but may not be physically accurate. A more detailed treatment requires X-ray binary population synthesis modeling for various IMFs (e.g., Fragos et al. 2013), which would involve variations in the mass ratio distribution with IMF. However, such treatments are likely to introduce many parameters into the analysis for which there are no plausible physical models for varying binarity, and for which there are no solid empirical constraints, making the benefit of such modeling inconclusive. Such work is beyond the scope of this paper (see Section 4).

Using Equation (8), and the constraints on ${N}_{38}/{L}_{K}$ presented in Figure 1(a), we computed maximum likelihood values for all 198,298 IMFs in our grid using the Cash statistic (cstat; Cash 1979). This procedure resulted in a normalized likelihood cube with three dimensions: ${\alpha }_{1}$, ${\alpha }_{2}$, and ${\gamma }_{{\rm{X}}}$. From our grid of models, the maximum likelihood model and 1σ errors are ${\alpha }_{1}={2.07}_{-1.1}^{+1.9}$, ${\alpha }_{2}={2.20}_{-0.24}^{+0.17}$, and ${\gamma }_{{\rm{X}}}={1.30}_{-0.45}^{+0.55}$ (10¹⁰ ${L}_{K,\odot }$)⁻¹ (see the 68% and 95% green contours in Figures 1(b)–(d)).

As we suspected, LMXBs provide stringent constraints on variations of the IMF for intermediate to massive stars (i.e., $m\gt 8$ ${M}_{\odot }$), and to the extent that we assume a single power law over the range of $m\gt 0.5$ ${M}_{\odot }$, a stringent constraint on ${\alpha }_{2}$. However, LMXBs essentially provide no useful constraints on the low-mass IMF slope.

Indeed, in Figure 1(a) we can surmise that the "invariant" case (red curve) is formally acceptable and within the 1σ threshold of our maximum probability model, while the interpolation to a single power-law IMF with $\alpha =2.8$ (blue curve) is not consistent with the LMXB observations. With an enhanced galaxy sample we not only confirm the result of P14 with better statistics, but also extend it by showing that it is possible to have an IMF that varies from being Kroupa for low-mass ellipticals to being bottom-heavy for high-mass ellipticals, as has been reported in the literature. In addition, the LMXB data suggest that the single power-law slope for such a varying IMF above 0.5 ${M}_{\odot }$ is unlikely to change much with velocity dispersion. In the next section, we utilize the formalism above, combined with a prior on how the IMF mass fraction due to low-mass stars can vary with velocity dispersion, to better constrain the values of ${\alpha }_{1}$, ${\alpha }_{2}$, and ${\gamma }_{{\rm{X}}}$.

As discussed in Section 3, LMXBs provide strong constraints on the variation of the IMF with velocity dispersion for stars with $m\gt 0.5$ ${M}_{\odot }$; however, variations at the low-mass end of the IMF are not well constrained by LMXBs alone. The green contours in Figures 1(b)–(d) show our 1σ and 2σ confidence intervals for ${\alpha }_{1}$, ${\alpha }_{2}$, and ${\gamma }_{{\rm{X}}}$ using only LMXB data.

It is clear that these data are unable to constrain well ${\alpha }_{1};$ however, we show that the model slope above 0.5 ${M}_{\odot }$, ${\alpha }_{2}$, and normalization factor, ${\gamma }_{{\rm{X}}}$, were strongly constrained. In order to better constrain ${\alpha }_{1}$, we utilized a prior derived by La Barbera et al. (2013), who concluded that the gravity-sensitive near-IR absorption features observed in the spectra of high-mass elliptical galaxies ($\sigma \ \approx$ 300 km s⁻¹) require ≈70%–90% of the total initial stellar mass to be contained in stars with $m\lt 0.5$ ${M}_{\odot }$, as integrated directly from the IMF. From Equation (1), we can compute the low-mass fraction for each of the n IMF models defined above following

$$\begin{eqnarray}{f}_{i}(\lt 0.5\,{M}_{\odot }) & = & \displaystyle \frac{M{(\lt 0.5{M}_{\odot })}_{i}}{{M}_{\odot }}\\ & = & \displaystyle \frac{{N}_{0,i}{2}^{{\alpha }_{2,i}-{\alpha }_{1,i}}}{{M}_{\odot }}{\displaystyle \int }_{0.1}^{0.5}{m}^{-{\alpha }_{1,i}+1}{dm}.\end{eqnarray} \tag{ 9 }$$

Using the low-mass fractions calculated via Equation (9), we next assigned a flat prior of 1 for $0.7\lt {f}_{i}\lt 0.9$, and 0 elsewhere for our grid of IMF models. Multiplying this prior by our likelihood cube (see Section 3), and renormalizing, resulted in a posterior probability distribution that we display as black contours in Figures 1(b)–(d). The resulting posterior probability provides a stringent constraint on ${\alpha }_{1}$ for massive ellipticals, resulting in best model values of ${\alpha }_{1}={3.84}_{-0.48}^{+0.09}$, ${\alpha }_{2}={2.14}_{-0.85}^{+0.20}$, and ${\gamma }_{{\rm{X}}}={1.30}_{-0.35}^{+0.65}$ (10¹⁰ ${L}_{K,\odot }$)⁻¹. The white curve, with gray 1σ envelope, displayed in Figure 1(a) shows our best model, which varies smoothly between a Kroupa IMF for low-mass ellipticals and a broken power-law IMF for high-mass ellipticals. The high-mass galaxy IMF component has a steep slope of ${\alpha }_{1}\,=$ 3.84 for stars $\lt 0.5$ ${M}_{\odot }$, and a slope of ${\alpha }_{2}\,=\,$ 2.14 for stars $\geqslant 0.5$ ${M}_{\odot }$. We note that while ${\alpha }_{2}$ is slightly flatter than the slope given by Kroupa, it still falls squarely within the uncertainty of the Kroupa IMF (${\alpha }_{2}=2.3\pm 0.7;$ Kroupa 2001).

By construction, this result is consistent with both the IMF for massive ellipticals being bottom-heavy, as inferred in the literature (see references in Section 1), and the IMF not varying significantly across the mass spectrum at the high-mass end, as inferred from the LMXB populations (see also P14). Additionally, the variable IMF has been constructed to yield a total (including stellar remnants) mass-to-light ratio that varies with velocity dispersion following the observational constraints from Cappellari et al. (2013). In absolute terms, the stellar mass-to-light ratio of our best-fit model varies from ${M}_{\star }/{L}_{K}=0.679$ at $\sigma =90$ km s⁻¹ to ${M}_{\star }/{L}_{K}=3.10$ at $\sigma =300$ km s⁻¹. The initial mass-to-light ratios (initial stellar mass over K-band luminosity) are $M/{L}_{K}=1.48$ at $\sigma =90$ km s⁻¹ to $M/{L}_{K}=4.12$ at $\sigma =300$ km s⁻¹, resulting in a ratio ${R}_{{\text{}}M/L}=2.78$. These values are below the limits placed on the dynamical values of mass-to-light ratio measured in the literature (e.g., Cappellari et al. 2013; Conroy et al. 2013).

We also noted in Section 1 that recent studies have found evidence for radial gradients in the IMFs of massive ellipticals, in which a bottom-heavy IMF appears to be appropriate for the inner regions of the galaxy, while a more Kroupa-like IMF is appropriate in the outer regions (see, e.g., Martín-Navarro et al. 2015; La Barbera et al. 2016a, 2016b). We tested to see whether the LMXB population showed evidence for such radial gradients by calculating the average ${N}_{38}/{L}_{K}$ values for LMXBs that were located within and outside the effective radii, r_e, of the population of massive elliptical galaxies. We utilized r_e values from Cappellari et al. (2013) to divide the ${L}_{X}\gt {10}^{38}$ erg s⁻¹ LMXB catalogs from Peacock et al. (2014) into inner- and outer-LMXB populations, and used the 2MASS K-band images to calculate the fraction of L_K that was within and outside r_e. We found consistent values of ${N}_{38}/{L}_{K}=1.4\pm 0.5$ and 1.7 ± 0.6 for the regions within and outside r_e, respectively, indicating that such gradients do not have a significant effect on our results.

To emphasize, our likelihoods are broadly insensitive to ${\alpha }_{1}$. We relied on La Barbera et al. (2013) as a meaningful prior and showed that we are able to reconcile their results with our analysis, with a statistically insignificant change in the derived values for the other two parameters, ${\alpha }_{2}$ and ${\gamma }_{{\rm{X}}}$. A choice of a different, meaningful prior would necessarily entail a different posterior probability density for ${\alpha }_{1}$; however, the essential constraint that this study places on ${\alpha }_{2}$ would not change significantly.

If a slight flattening of the high-mass slope of the generalized IMF is common in high-mass elliptical galaxies (shown to be required if we accept that the claims of a bottom-heavy IMF in the literature are correct), then several additional predictions result. We would expect a factor of ∼2 higher supernova rate in the high-redshift versions of these galaxies than their present-day stellar masses imply. A more severe effect might be seen on the rate of long-duration gamma-ray bursts (GRBs), which are often thought to come from only the most massive stellar explosions. At present, we have limited information about very high-redshift GRBs, and almost no information about very high-redshift supernovae (SNe). The information that exists suggests a GRB rate not easily explained by detectable star-forming galaxies (Tanvir et al. 2012), but the present results are easily explained by having most of the star formation take place in relatively small galaxies, below the HST detection threshold. We would expect similar increases in the rates at which double neutron star mergers produce short GRBs and gravitational wave sources, as well as an enhancement, albeit a somewhat lower one, of the rate of Type Ia SNe, since more white dwarfs would be produced and they would be skewed more heavily toward the massive end of the white dwarf spectrum.

It is also important to note that if the slope of the high-mass end of the IMF changes by about 0.1 dex, then the ratio of the number of core collapse supernovae from ≈40 ${M}_{\odot }$ stars to that from ≈8 ${M}_{\odot }$ stars changes by a factor of about 20%. The yields from the most massive core collapse supernovae can vary quite dramatically from those from the least massive ones (e.g., Kobayashi et al. 2006). As a result, abundances in the interstellar medium may provide a complementary test of the IMF. Because additional metal enrichment comes from thermonuclear supernovae, classical novae, and mass loss, a detailed treatment of this problem must be undertaken before a clear prediction can be made.

Future X-ray and high-resolution optical observations, and new population synthesis analyses could substantially improve constraints on the variation of IMF with velocity dispersion. In particular, new Chandra and HST observations of low-to-intermediate mass ellipticals would be helpful in ruling out a scenario where the high-mass end of the IMF is constant across all σ. Deeper Chandra observations of the low-mass galaxies in this study would have a similar effect, allowing us to probe to the more numerous population of low-luminosity LMXBs. In a forthcoming paper, M. Peacock et al. (2017, in preparation) will be examining the field population of LMXBs with ${L}_{{\rm{X}}}\gtrsim {10}^{37}$ erg s⁻¹ in low-mass elliptical NGC 7457 and comparing it with similar binaries detected in deep observations of massive ellipticals.

Testing the variations of the IMF in elliptical galaxies most effectively will require simultaneous modeling of near-IR spectroscopic data along with LMXB constraints using the combination of stellar population synthesis and X-ray binary population synthesis modeling (e.g., Fragos et al. 2013; Madau & Fragos 2016), in which the IMF, stellar ages, metallicities, and other physical parameters that influence X-ray binary formation are all modeled self-consistently. Such a population synthesis framework will be an important future step for advancing our knowledge of how the IMF varies among elliptical galaxies.

Additionally, a further understanding of the details of the stellar populations in the near-IR must be developed to ensure that the claims of a bottom-heavy IMF are reliable. At the present time, stellar evolution codes used for these purposes have limited or no treatment of unusual classes of stars such as interacting binaries and the products of binary interactions. Such stars are likely to be relatively unimportant for tracing out the optical bands where the models have been best calibrated. However, given that stars with the largest radius are most likely to interact, these stars may be increasingly important toward the reddest parts of the spectral energy distribution, where red giants and asymptotic giant branch stars are most important. For instance, it has already been shown that the S-type stars might appear with different frequencies in larger and smaller elliptical galaxies, and can potentially mimic the effects of bottom-heavy IMFs on the Wing–Ford band and the Na i D line, although not on Ca ii (Maccarone 2014).

Given the coincidence needed between changes at both the high- and low-mass ends of the IMF to reproduce the results we see here, as well as the subtlety of the features that have been used to suggest the bottom-heavy IMF, more work to investigate further the level of systematic effects from stars not included in standard stellar population synthesis models would be well motivated.

D.A.C., B.D.L., and R.T.E. gratefully acknowledge support from Chandra X-ray Center (CXC) grant GO4-15090A. D.A.C. acknowledges support from USRA and GSFC, as well as generous support from the University of Arkansas. A.K. acknowledges CXC grants GO4-15090C and GO5-16084B, and T.M. acknowledges GO4-15090B. M.P. and S.Z. acknowledge support from NASA ADAP grant NNX15AI71G, and M.P. additionally acknowledges CXC grant GO5-16084A and the Hubble Space Telescope (HST) grant HST-GO-13942.001-A.