MOSFIRE Spectroscopy of Quiescent Galaxies at 1.5 < z < 2.5. II. Star Formation Histories and Galaxy Quenching

We analyze the data using pyspecfit (Newman et al. 2014), a Bayesian code that fits stellar population templates to the observed spectroscopy and photometry of a galaxy. The code employs the nested sampling algorithm (Skilling 2004) to explore the multidimensional parameter space in an efficient way.

Our spectroscopic data typically cover the rest-frame wavelength range 3700–4200 Å, but the exact interval depends on redshift and target position on the MOSFIRE slit mask. Before the fit, in each spectrum, we mask out the pixels at the expected position of the [O ii] emission line, since our templates are purely stellar and do not include the contribution of ionized gas. To reduce the impact of sky lines on the spectrum, we also exclude from the fit the spectral pixels with the strongest sky emission, as measured by the observed sky spectrum.

We take the photometric data from the 3D-HST catalog, which contains publicly available measurements from the UV to the near-infrared. To avoid potential contamination from dust emission, we mask out the IRAC 8 μm channel. To account for systematic effects, we add in quadrature a 5% contribution to the uncertainty of each photometric point. Data points corresponding to nondetections are included in the fit with their formal uncertainty.

One limitation of optical and near-infrared spectra is the large uncertainty on the absolute flux calibration, which can be affected by slit loss and variable weather conditions. Instrumental effects can also change the overall shape of the spectral response, particularly near the limits of the observed wavelength range. Since the absolute calibration of photometric measurements is typically much more reliable, during the fit we change the overall shape of the observed spectrum, using a fourth-degree polynomial, so that it agrees with the shape given by the broadband photometry. This procedure does not affect any of the spectral features of interest, which extend over narrow spectral ranges. We also scale the error spectra to ensure that the spectroscopic and photometric data are given appropriate weighting in the fit. To do this, for each galaxy, we run a first fit with the sole purpose of calculating the χ², which we then use to rescale the error spectrum. The scaling factors are in the range 1.7–2.3, suggesting that the uncertainty on the spectra that is produced by the data reduction is indeed too small.

We generate the synthetic spectra using the Bruzual & Charlot (2003) library and adopting a Chabrier (2003) initial mass function (IMF) and the dust extinction law from Calzetti et al. (2000). A simple stellar population (SSP) of a given age can be used to generate a realistic spectrum once the following four parameters are specified: redshift, stellar velocity dispersion, dust attenuation A_V, and stellar metallicity Z. We leave these four parameters free during the fit, setting their priors in the following way. The spectroscopic redshift can be determined by eye with reasonable accuracy, and we use a narrow Gaussian prior centered on the values reported in Paper I. We also adopt a Gaussian prior for the stellar metallicity, centered on the solar value (Z = 0.02) but wide enough to encompass all reasonable values (σ_Z = 0.005). For the other parameters, we use uninformative priors in the form of uniform distributions: between 50 and 500 km s⁻¹ for the velocity dispersion and between 0 and 4 for A_V.

An additional free parameter is the age of the stellar population, which means that a model can be fully specified using five parameters. However, the assumption of SSP is likely too simplistic for real galaxies, and more complex SFHs, with a higher number of free parameters, are often considered. We attempt here a systematic exploration of different SFH models, with a twofold goal: first, we want to derive the stellar population properties in a robust way, ensuring that the results are not affected by the assumed form of the SFH; second, we are going to assess whether some of the SFH models are better than others in describing the spectra of real galaxies.

Each SFH model is simply a function of look-back time that is normalized to unity and depends on one or more parameters. The SFH models used in our analysis, together with the priors for each parameter, are listed in Table 1. We use the following six SFH models.

• 1.  
  The SSP represents the basic form of SFH, in which all stars form in an instantaneous burst. The only free parameter is the age of the population, for which we adopt a log-uniform prior that is truncated at the age of the universe at the observed redshift.
• 2.  
  Next, we consider a top-hat model, i.e., a burst with a finite duration. This model has two free parameters: the age and the duration of the burst.
• 3.  
  The most popular SFH choice is an exponentially declining function, also called a tau model. This SFH has two free parameters: the age and the timescale τ, and the star formation rate is proportional to $\exp (-t/\tau )$. In this SFH, the current star formation rate is lower than at any time in the past, which is a reasonable assumption for galaxies that are, by selection, quiescent.
• 4.  
  Another popular choice is the delayed tau model, which has the form $t\exp (-t/\tau )$. This SFH has been proposed in order to overcome the unrealistically sharp peak at early times featured by the tau model without introducing additional free parameters.
• 5.  
  In order to allow more flexibility, we introduce a new model that consists of a constant SFH followed by an exponential decay, hereafter called the constant + tau model. This model has three free parameters—age, duration of the flat phase, and timescale τ for the exponential phase—and can be considered a generalization of both the top-hat and the tau models. Its main advantage is that it allows the quenching phase to be completely detached from the initial star-forming phase.
• 6.  
  Finally, we consider a nonparametric SFH with a large degree of flexibility in order to explore further possibilities, such as secondary bursts or rejuvenation episodes. We adopt a piecewise description of the SFH with seven nodes at fixed, logarithmically spaced look-back times (0, 0.25, 0.45, 0.8, 1.4, 2.5, and 4 Gyr). Since the SFH models must be normalized, this corresponds to six free parameters. We choose to leave the amplitude of the first six nodes free, while the amplitude of the oldest node is arbitrarily fixed to unity. To avoid biasing the results, we keep the oldest node fixed at 4 Gyr for all galaxies, although this means that, according to the redshift, stars slightly older than the age of the universe may be allowed. The SFH is then constructed by linearly interpolating between the nodes. We adopt a lognormal prior for each amplitude, centered on 1 and with a width of 2 dex, thus allowing a large amount of freedom.

The SFH represents the relative contribution of stars of different ages to a population. The total luminosity is then determined by the total number of stars in the population, i.e., by the total stellar mass. For each template, this can be calculated by multiplying its mass-to-light ratio by the observed luminosity, corrected for cosmological dimming. The stellar mass is therefore an output of the fitting procedure but does not represent a free parameter of the model. The same is true for the instantaneous star formation rate, which can be easily derived from the SFH after accounting for the mass loss due to winds and supernovae. We make the simplifying assumption of instantaneous recycling, so that a fixed return fraction R of the newly formed stars is immediately lost and does not contribute to the stellar mass. For a Chabrier IMF, the return fraction is R = 0.41 (Madau & Dickinson 2014). The relation between the SFH and star formation rate at any cosmic time t is therefore

$$\begin{eqnarray}&&\mathrm{SFH}(t)\cdot {M}_{* }=\displaystyle \frac{{{dM}}_{* }(t)}{{dt}}=(1-R)\cdot \mathrm{SFR}(t).\end{eqnarray} \tag{ 1 }$$

Figure 1 illustrates the spectral fitting for one representative galaxy drawn from our sample. The top two panels show, in black, the observed photometry and the MOSFIRE spectrum. Six colored lines, corresponding to the best-fit templates for the six types of SFH, are overplotted. For each SFH model, the colored line in the top panel is exactly the same as the one in the second panel, because we simultaneously fit the spectroscopic and photometric data. In the bottom panel of the figure, the six types of SFHs are plotted in units of star formation rate, using Equation (1) for the conversion. In each panel, the best-fit SFH is shown as a colored line, while thin gray lines illustrate an even sampling of the posterior distribution. We also construct the median SFH by taking the median of the gray lines at each look-back time. This is shown as a thick black line and is more representative of the posterior distribution than the single best-fit curve. Throughout this work, when we quote the results of spectral fitting such as ages, star formation rates, etc., we always use the median calculated from an even sampling of the posterior distribution, not the best-fit value.

Zoom In Zoom Out Reset image size

We are able to obtain meaningful fits using each of the six SFH models for every galaxy in our sample except for one, ID 35616. Despite the fact that the MOSFIRE spectrum of this galaxy has a remarkably high signal-to-noise ratio, none of the six models of SFH is able to reproduce the photometry; for example, adopting a tau model, we obtain a reduced χ² (for the photometry only) of 4.3, which is a clear outlier given that the rest of the sample has a mean χ² of 1.1. Since the HST imaging shows that this galaxy is clearly interacting with a close companion (see Figure 2 of Paper I), it is likely that some of the photometric measurements are contaminated. Furthermore, this object is also the only one in the sample with strong Hα emission, as we discuss below in Section 3.5. For these two reasons, we exclude 35616 from the remainder of our analysis.

The 23 objects for which we obtain a good spectroscopic and photometric fit are shown in Figure 2. For each galaxy, we plot the MOSFIRE spectrum with the best-fit template, and the median SFH that we obtain when adopting a constant + tau model. The galaxies are sorted by their median age, from the youngest (top) to the oldest (bottom). The sample clearly spans a variety of stellar ages. A few galaxies, for which the fit indicates ages smaller than 1 Gyr, present strong higher-order Balmer lines, indicative of populations dominated by B and A stars. The vast majority of the galaxies in the sample, however, show moderate Balmer absorption lines and a relatively strong 4000 Å break, indicating the presence of both young and old stars. These galaxies can easily be sorted in age by looking at the relative strength of the two Ca ii absorption lines (see, e.g., Sánchez Almeida et al. 2012). If a significant population of stars younger than 1 Gyr is present, then the blend of the Ca H and H lines is deeper than Ca K. A further test of the relative ordering of these spectra comes from the Hδ absorption line, which is stronger for younger populations. This qualitative assessment shows that the results of the fits are consistent with the rest-frame optical spectra and not dominated by the photometry.

Zoom In Zoom Out Reset image size

In addition to the SFH, our spectral fits also determine the stellar mass, dust attenuation, and stellar metallicity, which are listed in Table 2 for the entire sample. The stellar metallicities span the range between half solar and slightly extrasolar, but their posterior distributions are broad, meaning that this quantity cannot be robustly measured by the spectral fits. We choose to include it as a free parameter anyway to avoid biasing the results for the other properties.

Table 2.  Stellar Population Properties of the MOSFIRE Sample

ID	Field	z	t50a	AVa	Za	log ${M}_{* }$/${M}_{\odot }$ a	log sSFRa,b	SFR(Hα)c
(3D-HST)			(109 yr)	(mag)			(yr−1)	(${M}_{\odot }$ yr−1)
30475d	UDS	1.63	0.30 ± 0.06	0.85 ± 0.10	0.013 ± 0.004	10.71 ± 0.02	−11.0 ± 1.0	<8.9
19958	COSMOS	1.72	0.48 ± 0.04	0.73 ± 0.06	0.012 ± 0.002	10.64 ± 0.01	−14.2 ± 3.0	⋯
5681	COSMOS	2.43	0.61 ± 0.07	0.64 ± 0.09	0.017 ± 0.003	10.91 ± 0.02	−12.9 ± 3.6	<14.9
16629	COSMOS	1.66	0.81 ± 0.07	0.44 ± 0.05	0.017 ± 0.003	10.59 ± 0.01	−17.0 ± 5.7	<7.2
22802	UDS	1.67	0.83 ± 0.1	0.62 ± 0.09	0.015 ± 0.004	11.01 ± 0.02	−13.2 ± 4.9	<6.8
13083	COSMOS	2.09	1.1 ± 0.1	0.20 ± 0.10	0.022 ± 0.003	10.95 ± 0.03	−20.0 ± 9.8	<14.5
1769	COSMOS	2.30	1.1 ± 0.2	0.43 ± 0.10	0.022 ± 0.005	11.16 ± 0.03	−11.1 ± 0.2	⋯
24891	UDS	1.60	1.1 ± 0.1	0.55 ± 0.10	0.015 ± 0.004	10.87 ± 0.02	−10.4 ± 0.09	<6.3
30737d	UDS	1.62	1.1 ± 0.1	0.51 ± 0.09	0.019 ± 0.004	11.23 ± 0.02	−10.6 ± 0.08	<15.1
17361	COSMOS	1.53	1.1 ± 0.09	0.37 ± 0.07	0.011 ± 0.002	10.74 ± 0.02	−10.9 ± 0.1	<1.7
25526	EGS	1.75	1.2 ± 0.2	0.34 ± 0.10	0.015 ± 0.005	10.69 ± 0.03	−10.9 ± 1.3	⋯
17926	EGS	1.57	1.2 ± 0.1	0.33 ± 0.09	0.012 ± 0.002	10.85 ± 0.05	−20.1 ± 10.6	<2.0
17255	COSMOS	1.74	1.3 ± 0.2	0.31 ± 0.06	0.017 ± 0.004	10.81 ± 0.03	−11.3 ± 0.4	⋯
29352	UDS	1.69	1.5 ± 0.3	0.23 ± 0.07	0.012 ± 0.003	10.82 ± 0.05	−12.6 ± 5.5	⋯
43367	UDS	1.62	1.5 ± 0.3	0.68 ± 0.10	0.017 ± 0.004	11.06 ± 0.04	−10.4 ± 0.1	<9.2
37529	UDS	1.67	1.5 ± 0.2	0.58 ± 0.07	0.014 ± 0.004	10.95 ± 0.04	−10.8 ± 0.08	<7.9
17641	COSMOS	1.53	1.5 ± 0.2	0.28 ± 0.05	0.013 ± 0.003	10.61 ± 0.04	−11.5 ± 1.8	<3.0
11494	COSMOS	2.09	1.6 ± 0.2	0.26 ± 0.06	0.022 ± 0.003	11.49 ± 0.03	−11.8 ± 0.09	<9.9
22719	EGS	1.58	1.8 ± 0.5	0.40 ± 0.09	0.012 ± 0.003	10.96 ± 0.08	−11.7 ± 7.4	<1.0
7884	COSMOS	2.11	2.0 ± 0.3	0.38 ± 0.06	0.013 ± 0.004	11.35 ± 0.02	−10.7 ± 0.08	<9.6
32707	UDS	1.65	2.1 ± 0.3	0.31 ± 0.08	0.014 ± 0.004	11.20 ± 0.04	−11.5 ± 0.10	<13.4
17364	COSMOS	1.53	3.0 ± 0.4	0.13 ± 0.07	0.015 ± 0.003	10.94 ± 0.02	−11.5 ± 0.07	<1.7
17089	COSMOS	1.53	3.4 ± 0.2	0.30 ± 0.06	0.013 ± 0.002	11.45 ± 0.02	−12.0 ± 0.07	<2.9

Note. More properties, including coordinates and morphological and dynamical measurements, are available in Table 2 of Paper I.

^aDerived from spectral fitting adopting a constant + tau model. The value and uncertainty are given by, respectively, the median and standard deviation of the posterior distribution. ^bAveraged over the past 100 Myr. ^c3σ upper limit on star formation rate, measured from the dust-corrected Hα flux. ^dDetected in the X-ray.

Download table as:  ASCIITypeset image

The example in Figure 1 shows that the six SFH models can produce different results when applied to the same object, in terms of both the best-fit template and the recovered physical properties, such as the stellar age. One of the goals of the present analysis is to compare different SFH models and assess the extent to which the choice of the model affects the results of the spectral fit.

Since each SFH model has a different number and type of parameters, we need to define nonparametric quantities in order to be able to compare the results among one another. Following Pacifici et al. (2016), we define the characteristic time t_X as the look-back time at which X% of the stellar mass was already formed: ${\int }_{{t}_{X}}^{\infty }\mathrm{SFH}({t}_{\mathrm{lb}}){{dt}}_{\mathrm{lb}}=X/100$, where t_lb is the look-back time. We start by considering t₅₀, i.e., the median stellar age. We calculate the value of t₅₀ using the six SFH models and show the results for the entire sample in the top panel of Figure 3. Since the tau model is the most widely used, we take it as a reference and plot, for each galaxy, the median ages obtained with the other models versus the median age obtained from the tau model. Clearly, the SSP model systematically underestimates the median age. This is a well-known issue, due to the fact that the SSP model effectively returns the stellar population age weighted by light instead of mass. Interestingly, the other five models yield consistent results, showing that the galaxies in our sample span an order of magnitude in stellar ages, from 300 Myr to 3 Gyr. The delayed tau and constant + tau models are in extremely good agreement with the tau model, while the top-hat and piecewise models present a slightly larger scatter. The piecewise model shows an abundance of points around t₅₀ ∼ 1.5 Gyr, which is likely an artificial effect due to the positioning of the nodes (see Figure 1).

Zoom In Zoom Out Reset image size

Next, we consider the quantity t₁₀ − t₉₀, which is the time it took to assemble the central 80% of the stellar mass of a galaxy, and can be considered as a proxy for the length of the star formation activity. The bottom panel of Figure 3 shows that different SFH models yield wildly different estimates for t₁₀ − t₉₀. In particular, the top-hat produces very long bursts, often extending all the way to the limit set by the age of the universe. The tau model, instead, produces very short bursts, which is not surprising given its simple parameterization. The piecewise model yields intermediate results with large uncertainties, which are likely more realistic given its high degree of flexibility. The SSP yields zero by definition and is not shown in the figure. We conclude that the length of the SFH is not as robust as the median age, and its measurement depends strongly on the assumed SFH model.

Another important property of quiescent galaxies is their current level of star formation activity. For each SFH model, we divide the star formation rate (given by Equation (1)) by the stellar mass, thus obtaining the specific star formation rate (sSFR). The SSP and top-hat models, by definition, have zero current star formation rate. The delayed tau and constant + tau models are very similar to the tau model; therefore, we only compare the piecewise and tau models in the top panel of Figure 4. For most of the galaxies, we are only able to obtain upper limits, indicating that the spectra and photometry are not able to constrain the current sSFR. But if we use the sSFR averaged over the previous 100 Myr, shown in the bottom panel, then we are able to obtain measurements with meaningful error bars for most of the galaxies. Moreover, the piecewise and tau models are in excellent agreement, confirming that the spectra are indeed able to detect the level of star formation in the recent but not immediate past. The values of sSFR we obtain range between 10⁻¹² and 10⁻¹⁰ yr⁻¹, corresponding to mass-doubling times longer than three times the age of the universe at z ∼ 2, consistent with the quiescent nature of the sample.

Zoom In Zoom Out Reset image size

The most robust measurements are therefore the median age t₅₀ and the sSFR averaged over the last 100 Myr, which are listed for each galaxy in Table 2. When looking at parameters that are not directly related to the SFH, different models give consistent results. The stellar mass is particularly robust, with a difference between the piecewise and tau models of only 0.02 ± 0.04 dex (mean and standard deviation, respectively). Discrepancies in dust extinction (−0.04 ± 0.10 mag) and stellar metallicity (7% ± 15%) are also small.

We now ask the question of whether the data favor any of the six models over the others. To answer this, we take advantage of the nested sampling algorithm employed by pyspecfit, which, in addition to the posterior distribution, also returns the model evidence. This quantity is the integral of the likelihood, weighted by the prior, over the full parameter space, and its value is proportional to the support of the data for a given model. In practice, model evidence is only relative: the ratio of the evidence of two models can be interpreted as the ratio of the odds in favor of one model against the other. Once again, we take the tau model as reference, and for each galaxy, we define the quantity ${\rm{\Delta }}\mathrm{ln}\,\mathrm{evidence}$ as the logarithmic difference between the evidence for a given model and the evidence for the tau model.

We show the distribution of ${\rm{\Delta }}\mathrm{ln}\,\mathrm{evidence}$ for the entire sample as a function of galaxy age in Figure 5. This quantity can be directly translated to the odds that the tau model explains the data better than another model, and these are shown on the right-hand y axis. We also mark the regions of the plot corresponding to weak, moderate, and strong evidence according to a popular, though somewhat arbitrary, scale (see, e.g., Trotta 2008). We find that the delayed tau and constant + tau models fare as well as the tau model. On the other hand, the SSP and top-hat models are strongly disfavored. It is particularly interesting to note that the SSP actually fares slightly better than the tau model for young galaxies but is unable to well describe older galaxies, for which a more extended SFH is clearly needed. Similar conclusions apply to the top-hat model.

Zoom In Zoom Out Reset image size

The piecewise model is comparable to the tau model except at young ages, where it performs moderately worse, presumably because of its rather coarse binning. However, since the evidence is particularly sensitive to the adopted priors, the comparison of nonparametric and parametric models must be performed with care. For example, the maximum allowed stellar age for the piecewise model is fixed to 4 Gyr (see Section 3.1), while for the parametric models, we adopt the age of the universe at the observed redshift. To test whether this has any effect on the results, we run an additional spectral fit using a piecewise SFH where the last node, instead of being fixed at 4 Gyr, is equal to the age of the universe at the observed redshift. The resulting evidence is, on average, unchanged, and individual points move in Figure 5 by $\left|{\rm{\Delta }}\,\mathrm{ln}\,\mathrm{evidence}\right|\lt 1.5$. We also test for the width of the priors, which we vary by one order of magnitude in each direction (i.e., adopting LogNormal(0, 1) and LogNormal(0, 3) on each of the SFH nodes; see Table 1). The quantity ${\rm{\Delta }}\,\mathrm{ln}\,\mathrm{evidence}$ tends to decrease in both cases, on average by 1.0 and 0.5, respectively, and the result shown in Figure 5 is qualitatively unchanged. The physical properties, such as median ages and stellar masses, are remarkably stable. For an in-depth discussion of the effects of the prior when using nonparametric models, we refer to Leja et al. (2018).

Keeping in mind that the evidence is weighted by the size of the parameter space, we conclude that the extra degrees of freedom of the piecewise models (six free parameters, as opposed to the two of the tau model) are not warranted, since the final results are not measurably better than those obtained with simpler models. This means that the secondary bursts recovered by the piecewise model in some cases, such as the one shown in Figure 1, are possible but not required by the data. We caution, however, against a simplistic interpretation of this result: in the context of SFH models, having more degrees of freedom than required by the data is not necessarily a poor choice (see, e.g., Simha et al. 2014). For example, if one is interested in the measurement of a specific parameter, such as the median stellar age, it may be a good idea to allow more freedom than strictly required by the data.

Our analysis can be summarized as follows: the tau, delayed tau, and constant + tau models yield consistent results, whereas the SSP and top-hat models show clear limitations in reproducing the data. The piecewise model is the most flexible one but does not describe the data significantly better than the simple tau model. It seems that the constant + tau model strikes a balance between simplicity and flexibility, and for this reason, we choose it as our "fiducial" model.

Our additional, shallower H- and K-band spectroscopic observations target the Hα emission for 19 out of 24 galaxies in the sample. We find only one clear detection of Hα, for 35616, the object we removed from the main analysis. The [N ii] line is also detected, and, although Hα is contaminated by a strong sky line, we are able to derive a flux ratio [N ii]/Hα ∼ 0.8, which is high and suggests a contribution from active galactic nucleus (AGN) activity (e.g., Kewley et al. 2013), consistent with the fact that this is one of the three galaxies in the sample that are detected in the X-ray. About 10% of quiescent galaxies at high redshift harbor this type of emission (Belli et al. 2017a). Since 35616 is in a close interaction, it is possible that the AGN activity is triggered by gas that has been tidally stripped from the companion.

For the remaining galaxies in the sample, we scale the observed spectrum so that the median level of the continuum around Hα matches the best-fit template from the spectral fit to ensure a proper flux calibration. We then subtract the best-fit template to remove the contribution of the stellar absorption. We calculate upper limits on the nondetections by integrating the error spectrum in a window of 500 km s⁻¹ (corresponding to a line width σ ∼ 200 km s⁻¹) centered on the Hα wavelength. After correcting the line fluxes for the extinction due to dust, adopting the attenuation A_V obtained from the spectral fit, we convert to values of star formation rate (based on a Chabrier IMF) following Kennicutt (1998). The 3σ upper limits on the star formation rates are listed in Table 2. The upper limits range from 1 to 15 ${M}_{\odot }$ yr⁻¹; the highest values correspond to cases in which strong sky lines enter the Hα spectral window. For each galaxy, the limit on the star formation rate derived from Hα is higher than the value predicted by the spectral fit.

To obtain a deeper measurement, we shift each of the 18 spectra to the rest frame and stack them. The result is shown in the top panel of Figure 6: both Hα and [N ii] are clearly detected. The bottom panel shows what we obtain when stacking the spectra without subtracting the stellar template. In this case, the emission from [N ii] is detected, while Hα is marginally detected in absorption. Without a detailed template of the stellar continuum, it would not be possible to measure such a small amount of Hα emission.

Zoom In Zoom Out Reset image size

We model the emission lines Hα, [N ii] 6548, and [N ii] 6584 as three Gaussians. We fix the wavelengths to their theoretical values and the flux ratio of the two [N ii] lines to 3 and set the velocity dispersion of the three lines to be identical. The best fit of this model to the continuum-subtracted stack is shown in red in Figure 6. We obtain an Hα flux of (6.7 ± 0.8) · 10⁻¹⁸ erg s⁻¹ cm⁻², a flux ratio of [N ii]/Hα = 1.0 ± 0.1, and a velocity dispersion of 233 ± 17 km s⁻¹. One galaxy, 30475, has a marginally detected [N ii] line (it is also one of the three galaxies of the sample detected in the X-ray). Removing it from the stack has a minimal impact on the results. Further splitting the sample into two random subsets does not change the results of the fit in a statistically significant way, confirming that the ionized emission is not due to any individual object. We also test the effect of the SFH model on the subtraction of the stellar template: adopting the best-fit piecewise model instead of the constant + tau model, we obtain nearly identical emission-line measurements.

The high [N ii]/Hα ratio we obtain for the stacked spectrum rules out star formation activity as the main source of emission lines. In the local universe, quiescent galaxies commonly feature weak emission lines due to the ionization from old stellar populations (Cid Fernandes et al. 2011). However, the maximum Hα equivalent width predicted for this mechanism is 3 Å, while we measure 4.7 ± 0.6 Å in our stack. Therefore, we conclude that the most likely source for the ionized gas emission is low-level AGN activity. Ongoing star formation could still contribute to the observed flux, which therefore can be used to derive a strict upper limit to the star formation rate. After correcting the observed Hα flux for the average value of dust extinction in the sample (A_V = 0.43), we obtain an star formation rate upper limit of 0.9 ± 0.1 ${M}_{\odot }$ yr⁻¹, which confirms the truly quiescent nature of the sample. Dividing by the average stellar mass, we obtain $\mathrm{log}\,\mathrm{sSFR}\sim -11.1$, which is fully consistent with the values derived from the spectral fits and shown in Figure 4.

Our analysis of the rest-frame colors and number densities demonstrates that not all quiescent galaxies go through the post-starburst phase after being quenched. This implies a variety of quenching timescales: some galaxies undergo a fast quenching, becoming post-starburst and then evolving along the red sequence, while the remaining ones must follow a slow-quenching path that never produces the blue colors typical of post-starburst galaxies. We illustrate these two paths on the UVJ diagram in Figure 12 using two simple toy models. The fast-quenching track, shown in blue, is exemplified by a tau model with a short timescale, τ = 100 Myr, while for the slow track, in green, we take τ = 1 Gyr. We also include the effect of dust reddening, which we model in an identical way for the two tracks: the dust attenuation is A_V = 2 at the beginning, declines following the evolution of the star formation rate, and then reaches a plateau at A_V = 0.4, the typical value found for our quiescent galaxies.

Zoom In Zoom Out Reset image size

The toy models are effective at illustrating the salient features of the fast and slow tracks. Clearly, the red sequence on the UVJ diagram can be joined at any point along its diagonal edge, a region often called the green valley. The exact position of the crossing point will depend on both the quenching timescale and the dust attenuation. The qualitative difference between the two tracks is, however, mainly due to the quenching timescale. For τ = 1 Gyr, for example, even a dust-free model will not enter the post-starburst region. In order to reproduce the colors of the most extreme post-starburst galaxies, which are even bluer than our selection box, it is necessary to have little dust and extremely short timescales, τ < 100 Myr. The rarity of these very blue quiescent galaxies, evident in Figure 12, implies that galaxies move across that region very rapidly, consistent with an extrapolation of our age–color relation.

Our number density analysis leads to a quantitative result: in order to explain the observed color distribution of quiescent galaxies, the two quenching channels must be in place simultaneously and have similar quenching rates at z ∼ 2. This is consistent with the qualitative result of our earlier study of the LRIS sample (Belli et al. 2015), in which we proposed that post-starburst galaxies follow a fast- quenching path based on the derived SFHs and the observation that they have systematically smaller sizes than older quiescent galaxies. Other studies, at both low and high redshift, have reached similar conclusions (e.g., Schawinski et al. 2014; Moutard et al. 2016; Wild et al. 2016; Carnall et al. 2018; Forrest et al. 2018; Wu et al. 2018a). In these works, the two quenching channels are often characterized not only by different timescales but also by different morphological or structural properties, as we discuss below.

Given their importance in the context of galaxy quenching, we want to explore in more detail the properties of post-starburst systems. For this analysis, we make use of the catalog from the 3D-HST survey (Brammer et al. 2012), which covers a much smaller area than UltraVISTA but is based on high spatial resolution imaging and grism spectroscopy from HST. In order to explore a variety of physical properties, we assemble different catalogs that have been publicly released and include grism redshifts (Momcheva et al. 2016), rest-frame colors and stellar population properties (Skelton et al. 2014), size and morphological properties (van der Wel et al. 2014), and measurements of the local overdensities (Fossati et al. 2017). Given the limited survey area, the number of post-starburst galaxies is small; to increase the statistics, we slightly extend the mass range, selecting all objects with $\mathrm{log}{M}_{\odot }/{M}_{* }\gt 10.6$. Still, we find only four galaxies with rest-frame colors in the post-starburst region in the redshift range 1 < z < 1.5; therefore, we limit the analysis to the range 1.5 < z < 2.5, where the cosmic volume probed is larger. With this selection, we obtain a sample of 554 quiescent galaxies, of which 65 fulfill the post-starburst color selection.

We plot some of the physical properties for this sample of quiescent galaxies in Figure 13 as a function of the rest-frame color ${S}_{{\rm{Q}}}$, defined in Equation (2). Using the age–color relation, we can then convert this into an inferred stellar age, shown on the top axis. The post-starburst galaxies are those with an inferred age between 300 and 800 Myr. The figure shows that this population has clearly distinct properties compared to the older quiescent systems. Post-starburst systems tend to have lower stellar masses, consistent with the mass function study performed by Wild et al. (2016). To compare the structural properties, we calculate a mass-normalized size for each object by taking the measured effective radius and dividing it by the size typical for quiescent galaxies at that stellar mass and redshift, which we take from the van der Wel et al. (2014) mass–size relation. The distribution of normalized sizes shows that post-starburst galaxies are more compact than the overall quiescent population. This was first shown by Whitaker et al. (2012) and successively confirmed by a number of studies (Belli et al. 2015; Yano et al. 2016; Almaini et al. 2017; Maltby et al. 2018; Wu et al. 2018a). The distribution of Sérsic indices does not vary strongly with stellar age, with a marginal increase toward the post-starburst region, indicating that these systems are likely spheroidal (Almaini et al. 2017; Maltby et al. 2018). Finally, in the last panel, we show the distribution of the local overdensity, finding that roughly 80% of post-starburst galaxies are found in an overdense environment (δ > 0), as opposed to 25% of the older quiescent galaxies. Despite residing in overdensities, post-starburst galaxies are not likely to be satellites: using the Fossati et al. (2017) catalog, we find that three quarters of the post-starburst galaxies shown in Figure 13 are centrals with a probability larger than 90%. Our results are consistent with the excess of post-starburst galaxies that has been found in clusters at z < 1 (Moutard et al. 2018; Socolovsky et al. 2018). However, these post-starburst galaxies at lower redshifts are less massive ($\mathrm{log}{M}_{* }/{M}_{\odot }\lt 10.5$) and are likely to be satellites.

Zoom In Zoom Out Reset image size

The emerging picture is one in which the fast-quenching channel produces post-starburst galaxies that are preferentially compact and spheroidal, have lower stellar masses, and lie in higher density environment compared to the bulk of the quiescent population. The fact that post-starburst galaxies have such distinct physical properties points to a physical difference between fast and slow quenching, rather than to a single physical process characterized by a wide distribution of quenching timescales. One pressing question remains: What are the physical mechanisms responsible for the two quenching channels?

Prime candidates for the fast-quenching mechanism are gas-rich processes that naturally explain both the short timescales and the compact sizes found in post-starburst galaxies. Major mergers (Hopkins et al. 2006) or violent disk instabilities (Dekel & Burkert 2014) can drive nuclear inflows and cause centrally concentrated starbursts. These processes have been explored in detail in recent theoretical work, specifically because they can explain why massive quiescent galaxies at high redshift are generally very small compared to their local counterparts. Using a suite of hydrodynamical simulations, Zolotov et al. (2015) showed that mergers or violent disk instabilities can trigger a compaction event, in which a massive stellar core is formed in a short time, with a resulting decrease in the galaxy effective size. Quenching follows naturally due to a combination of gas exhaustion and stellar feedback; AGN feedback may also be involved, particularly for the maintenance of quiescence over longer timescales. However, these simulations predict a combined timescale for compaction and quenching of the order of 2 Gyr (see also Tacchella et al. 2016), which is substantially longer than what we measure for the post-starburst objects and more similar to our slow-quenching channel. This discrepancy may indicate the need for stronger feedback or additional physical mechanisms. A hint about the important role played by AGN feedback in this critical phase comes from the high incidence of ionized gas emission (not due to star formation) found in the post-starburst region of the UVJ diagram (Belli et al. 2017a).

Compaction can be triggered by a number of different processes, including major or minor mergers, counterrotating streams, perturbation due to giant clumps, and tidal interactions (Zolotov et al. 2015). The large fraction of post-starburst galaxies residing in overdensities suggests a prominent role for mergers and interactions with nearby galaxies (see also, e.g., Brodwin et al. 2013; Wu et al. 2014). We note that post-starburst systems are mostly centrals and therefore not affected by so-called environmental quenching, which has been proposed for explaining the abundance of quiescent low-mass systems in dense environments, and applies only to satellites (Peng et al. 2010). A link with environmental quenching, however, cannot be ruled out at z < 1, where the post-starburst population consists of lower-mass galaxies (Maltby et al. 2018).

Compaction is not the only possible explanation for the small sizes of quiescent galaxies at high redshift. The alternative possibility is that these systems formed at very high redshift, when star-forming disks were significantly smaller, and then evolved passively retaining their compact structures (Wellons et al. 2015). However, the young ages we measure for post-starburst galaxies rule out this scenario. We conclude that the observed properties of post-starburst galaxies require a compaction phenomenon, defined in its broadest sense as the formation of quiescent galaxies that are more compact than their star-forming progenitors. To establish exactly how and why this happens, more theoretical and observational studies are clearly needed.

In Section 5, we found that the importance of the fast-quenching channel increases with redshift, up to z ∼ 2.5. At even higher redshift, it is likely that fast quenching represents, in practice, the only available track, simply because of the young age of the universe. The earliest quiescent galaxy with a spectroscopic confirmation is a massive system at z = 3.717 that features strong Balmer absorption lines, typical of post-starburst systems (Glazebrook et al. 2017). Interestingly, this object happens to have a massive close companion (Simpson et al. 2017; Schreiber et al. 2018b). It is possible that the very first massive quiescent galaxies, which reside in overdense environments, were all quenched rapidly via major mergers.

Normal star-forming galaxies have depletion times of about 1 Gyr (e.g., Genzel et al. 2010), which means that a continuous replenishment of the gas reservoir via gas inflows is needed to continue the star formation activity over a cosmological timescale. If such replenishment is halted for any reason, the galaxy slowly consumes the remaining cold gas and quenches over long timescales. Such gas starvation is therefore the simplest explanation for the slow-quenching path. Unfortunately, there are many different physical mechanisms that can produce gas starvation. The most commonly invoked scenario is one where the gas cannot cool when the halo grows beyond a threshold mass due to the formation of a stable virial shock (Birnboim & Dekel 2003; Kereš et al. 2005), possibly in conjunction with heating by AGN radio-mode feedback (Bower et al. 2006; Croton et al. 2006). There are, however, other possibilities: the gas may cool but then be made stable against collapse by the deep potential well of a large stellar bulge (Martig et al. 2009; Genzel et al. 2014), or the gas inflows could stop even before reaching the halo (Feldmann & Mayer 2015).

In the local universe, the stellar metallicity of quiescent galaxies is substantially larger than that of star-forming objects of the same mass, which is a signature of starvation (Peng et al. 2015). This is consistent with the increasing importance of the slow-quenching channel at later cosmic epochs.

Another way to constrain the quenching mechanisms is to link quiescent galaxies with their star-forming progenitors. Using the SFHs from our spectral fits, we can estimate physical properties such as stellar mass and star formation rate for the progenitors of the quiescent galaxies in our sample. We caution that this analysis is strongly affected by the assumption of the SFH model and should be considered as purely qualitative.

In Figure 14, we show the inferred past evolution in both star formation rate (top) and stellar mass (bottom) for our MOSFIRE sample, split into two redshift bins. The tracks shown in the figure correspond to the constant + tau model for the SFH and are color-coded according to the median stellar age. We compare the inferred evolution with the properties of spectroscopically confirmed submillimeter galaxies (Toft et al. 2014) and compact star-forming galaxies (van Dokkum et al. 2015), two populations that have been proposed as possible progenitors of massive quiescent galaxies. The masses and star formation rates of compact star-forming galaxies are consistent with those of the progenitors of quiescent galaxies, at least for those observed at 1.5 < z < 2, and are also consistent with the star-forming main sequence (e.g., Whitaker et al. 2014). On the other hand, the submillimeter galaxies have substantially higher values of star formation rate, which are not consistent with the SFH tracks. However, this discrepancy can be easily explained by the typical duty cycle of submillimeter galaxies, which may be as short as 40 Myr (Toft et al. 2014); such a short-lived phase would leave virtually no signature on the spectra observed more than a Gyr later. Figure 14 also shows that the peak value of star formation rate is roughly the same for all quiescent galaxies, around 100–300 ${M}_{\odot }$ yr⁻¹, and does not depend on the median stellar age. This means that the SFHs of post-starburst galaxies do not require a burst of star formation above the main sequence; they simply require rapid formation and quenching.

Zoom In Zoom Out Reset image size

In the compaction scenario, galaxies first shrink and then quench. The population of compact star-forming galaxies is thought to represent the intermediate stage after compaction and before quenching (Barro et al. 2013; van Dokkum et al. 2015). Measurements of the kinematics have been used to draw a link between compact star-forming galaxies and quiescent systems. The integrated velocity dispersion of the gas in compact star-forming galaxies is in broad agreement with that of the stellar populations in quiescent galaxies (Barro et al. 2014; Nelson et al. 2014). A comparison of the spatially resolved kinematics, though, would be even more compelling. Recent observations showed that the gas in compact star-forming galaxies is actually rotating rapidly in both the ionized (van Dokkum et al. 2015; Wisnioski et al. 2018) and molecular (Barro et al. 2017; Tadaki et al. 2017) components. This is not in agreement with the kinematics of the most massive quiescent galaxies at z ∼ 0, which tend to be pressure-supported. At high redshift, however, quiescent galaxies seem to be very different from local ellipticals.

Our previous analysis of the MOSFIRE sample, which focused on the stellar kinematics, suggested that the rotational support of quiescent galaxies is significantly larger at high redshift (Paper I). This surprising result is consistent with detailed studies of a few gravitationally lensed systems at z > 2 (Newman et al. 2015, 2018b; Toft et al. 2017) and a study of a large sample at z < 1 (Bezanson et al. 2018). The fact that massive quiescent galaxies at high redshift present a large degree of rotational support has two important consequences in the context of quenching. First, it suggests that the giant ellipticals observed at z ∼ 0 were originally quiescent disks, which lost their angular momentum likely because of dry mergers. This means that, contrary to the standard view, the morphological transformation does not need to happen simultaneously with the quenching process in order to form ellipticals. Second, it offers an additional constraint on the quenching process, which needs to produce a rotating stellar system. Since starvation leaves the morphology and kinematics untouched, the slow-quenching channel should produce preferentially rotating quiescent systems. However, gas-rich mergers can also produce rotating remnants (Wuyts et al. 2010). It is therefore not clear whether kinematics can be used for a clean identification of the progenitors.

The question of whether the two quenching channels require distinct populations of progenitors is still open. It is possible that all star-forming galaxies can potentially undergo fast quenching after a close interaction or merger. On the other hand, post-starburst galaxies show smaller stellar masses compared to older quiescent systems, and the tracks shown in Figure 12 also indicate a mass segregation between fast and slow quenching. If this were the case, star-forming galaxies above a certain stellar mass (which depends on redshift) would be the ideal candidates to be progenitors of the slow-quenching track.

We analyzed the deep near-infrared spectra of 24 massive quiescent galaxies at 1.5 < z < 2.5 obtained with Keck/MOSFIRE. The high quality of the spectroscopic data, which we now release publicly, allowed us to perform a study of the stellar population properties of these objects to an unprecedented level of detail.

We took full advantage of all of the available data by fitting templates simultaneously to the observed spectroscopy and broadband photometry. We explored six different SFH models with the goal of understanding how this impacts the derived stellar population properties. We found that the SSP and top-hat models are not able to reproduce the observed data for at least some galaxies in our sample. The remaining four models were all equally successful, including the commonly used tau model, which has only two free parameters. Despite the degeneracy among different SFH models, the median stellar ages are found to be robust, together with the sSFR averaged over the past 100 Myr. However, the duration of the star formation depends strongly on the SFH shape and therefore cannot currently be measured robustly.

In order to break the degeneracy between different SFH models, significantly deeper spectra would be needed. This may be possible for individual objects, but for larger samples, we will have to wait for the next generation of ground-based telescopes. A complementary approach is to use the stellar metallicity to indirectly probe the formation timescale (Kriek et al. 2016), which, however, requires a wider wavelength range than that available for our sample.

By stacking shallower spectra obtained at redder wavelengths, we detected the Hα and [N ii] emission lines. The Hα line is faint and would not be detected without a proper subtraction of the underlying stellar absorption. We measured a high line flux ratio, [N ii]/Hα = 1.0 ± 0.1, which rules out star formation as the sole origin of the ionization. About 10% of quiescent galaxies at high redshift feature emission lines with a similar flux ratio, which is likely due to AGN-driven outflows (Belli et al. 2017a). Our finding suggests that low-level AGN activity is even more common in quiescent galaxies but in most cases is too weak to be noticed. This is also confirmed by the sample of gravitationally lensed quiescent galaxies in Newman et al. (2018a), where the lensing magnification allows us to observe faint emission lines in most of the objects. Such low-level AGN feedback may play a role in either quenching the star formation or maintaining the galaxies as quiescent.

Assuming that the dust-corrected Hα flux is entirely due to star formation yields a strict upper limit on the star formation rate, which is 0.9 ± 0.1 ${M}_{\odot }$ yr⁻¹. This corresponds to an average $\mathrm{log}\,\mathrm{sSFR}\sim -11.1$, confirming the results of the spectral fitting. This low value of star formation is in broad agreement with other studies that used less-direct methods (Fumagalli et al. 2014; Man et al. 2016a; Gobat et al. 2017) and suggests a very effective quenching in massive galaxies already at z ∼ 2.

By combining the MOSFIRE observations with our previous LRIS sample at lower redshift, we found that the majority of quiescent galaxies at 1 < z < 2.5 lie on a remarkably tight relation between stellar age and UVJ colors, which is independent of redshift and stellar mass. Such a relation deviates from that predicted by a simple model of passive evolution due to the contribution of dust reddening, which evolves rapidly during the post-starburst phase. The use of a large sample of high-quality spectra is therefore critical for determining the age–color relation. We quantified this relation and its scatter, which is only 0.13 dex in age, and we gave an explicit calibration that can be used to infer the approximate stellar age of any high-redshift quiescent galaxy, given its rest-frame colors, without the need for spectroscopic data. This is a simple yet powerful method that allowed us to model the evolution of the quiescent population and that can also be used for other applications, such as selecting specific targets for spectroscopic follow-up.

Based on our observations and on previous studies (e.g., Whitaker et al. 2015; Merlin et al. 2018), we argued that the lower and upper edges commonly used to define the quiescent region in the UVJ diagram are not justified and artificially exclude very young and very dusty quiescent galaxies. The young quiescent galaxies, also called post-starburst, represent a key evolutionary phase in the formation of the red sequence. The nature of the dusty objects (which are the only galaxies in our sample that do not lie on the age–color relation), instead, is currently not known: it is possible that they consist of edge-on disks and/or interacting systems.

It would be interesting to extend the age–color relation to higher redshifts; however, the few spectroscopic measurements available at z > 3 show that UVJ colors may become less effective in dividing quiescent and star-forming galaxies (Schreiber et al. 2018a).

Combining the age–color relation, which we derived from our spectroscopic sample, with a large, volume-limited sample drawn from the UltraVISTA survey, we were able to model the evolution of the red sequence at 1 < z < 2.5 in a quantitative and self-consistent way. The main result of this analysis is the identification of two quenching channels: fast quenching produces compact post-starburst objects, while slow quenching contributes to the growth of the red sequence by adding objects with stellar populations that are already old. Fast quenching is increasingly more important at high redshift, contributing about half of the growth of the quiescent population at z ∼ 2.

Post-starburst objects are clearly distinct from the bulk of the red sequence and feature lower stellar masses, compact sizes, spheroidal morphologies, and higher overdensities, as also shown by a number of recent studies. Together with their rapid formation timescale, these properties well fit a scenario in which violent, gas-rich events such as mergers are responsible for their formation. The slow-quenching channel, instead, must be due to a smooth decline in the gas reservoir, which may be caused by fuel starvation. However, there are many physical processes that may be responsible for starvation, including AGN feedback, virial shock heating, and cosmological or gravitational effects.