Stellar Mass Growth of Brightest Cluster Galaxy Progenitors in COSMOS Since z ~ 3

TL;DR

This study investigates the evolution of brightest cluster galaxy progenitors from redshift 3 to present, highlighting the changing roles of star formation and mergers in their mass growth across cosmic time.

Contribution

It provides a detailed analysis of the stellar mass growth phases of BCG progenitors, emphasizing the transition from in situ star formation to dry mergers and environmental effects.

Findings

BCG progenitors undergo three evolutionary stages.

Star formation dominates at high redshift ($z > 2.25$).

Dry mergers are key at low redshift.

Abstract

We examine the role of environment on the in situ star formation (SF) hosted by the progenitors of the most massive galaxies in the present-day universe, the brightest cluster galaxies (BCGs), from $z \sim 3$ to present in the COSMOS field. Progenitors are selected from the COSMOS field using a stellar mass cut motivated by the evolving cumulative comoving number density of progenitors within the Illustris simulation, as well as the Millennium-II simulation and a constant comoving number density method for comparison. We characterize each progenitor using far-ultraviolet--far-infrared observations taken from the COSMOS field and fitting stellar, dust, and active galactic nucleus components to their spectral energy distributions. Additionally, we compare the SF rates of our progenitor sample to the local density maps of the COSMOS field to identify the effects of environment. We find that BCG progenitors evolve in three stages, starting with an in situ SF dominated phase ($z > 2.25$). This is followed by a phase until $z \sim 1.25$ where mass growth is driven by in situ SF and stellar mass deposited by mergers (both gas rich and poor) on the same order of magnitude independent of local environment. Finally, at low redshift dry mergers are the dominant stellar mass generation process. We also identify this final transition period as the time when progenitors quench, exhibiting quiescent NUV\emph{rJ} colors.

Tables

Table 1. Table 1 Illustris-selected BCG Progenitor Sample Distribution

$z$	Total Progenitors	Total Photo-$z$	Total Spec-$z$	Total AGN†	X-ray AGN	MIR AGN	Radio AGN	SED3FIT AGN‡
0.5-0.8	26	2	24	2 (7.7%)	1	1	2	0
0.8-1.1	180	74	106	14 (7.8%)	1 2	6	0	5
1.1-1.5	216	156	60	8 (3.7%)	4	4	1	3
1.5-2.0	401	371	30	29 (7.2%)	22	15	1	12
2.0-2.5	347	308	41	32 (9.2%)	15	24	1	9
2.5-3.0	274	254	20	31 (11.3%)	10	26	0	8

Equations

$$⟨N(z){⟩}^{\prime }=N_0+\Delta z(A_0^{\prime }+A_1^{\prime }{N}_0+A_2^{\prime }{N}_0^2)$$

$$+\Delta z^2(B_0^{\prime }+B_1^{\prime }{N}_0+B_2^{\prime }{N}_0^2),$$

Full text

††slugcomment: Received 2019 May 8; revised 2019 June 29; accepted 2019 July 8; published 2019 August 22

1. Introduction

The evolutionary history of today’s most massive galaxies is an important component in understanding the evolution of large-scale structure in the universe and galaxy populations in general, placing vital limits on the speed that baryonic matter can assemble into gravitationally bound structures. These are massive elliptical galaxies, called brightest cluster galaxies (BCGs), and are normally located at the center of their cluster in position/velocity space (Lauer et al., 2014). They are composed of an old stellar population and a hot gaseous component in an extended envelope larger than predicted by the distribution of neighboring non-BCG ellipticals (Oemler, 1976). While most are quiescent, those that inhabit relaxed, massive haloes may exhibit atypically large ($>1\;M_{\odot}$ yr*-1*) star formation rates (SFRs). This is due to the surrounding intra-cluster medium (ICM) cooling and precipitating into the BCGs, forming cool-core galaxy clusters (e.g., Fabian, 1994; O’Dea et al., 2008, 2010; Rawle et al., 2012; Fogarty et al., 2019), or due to the accretion of gas-rich satellites (e.g., Tremblay et al., 2014; Webb et al., 2015a; Cooke et al., 2016). BCG stellar mass and SFRs may also be affected by the mass of their cluster halos at low redshift (e.g., Cooke et al., 2018; Gozaliasl et al., 2018). Their atypically high mass, quiescence, and extended light profiles collectively point to a formation mechanism unique with regard to their satellites, in which most of their stars formed at $z>2$, and conglomerated into a single BCG through many major and minor mergers (De Lucia & Blaizot, 2007). With this complex model in mind, a full identification of BCG progenitors requires the characterization of all galaxies at high redshift that are likely to merge and contribute their stellar mass onto proto-BCG “seeds”.

BCGs inhabit the high-mass end of the galaxy stellar mass function (SMF); at $z\sim 0.5$, BCGs in massive, relaxed clusters are typically $\sim 10^{11.5}\;M_{\odot}$ or greater (e.g., Lidman et al., 2012; Lin et al., 2013; Chiu et al., 2016; Cooke et al., 2016). At this stellar mass, the population is dominated by passive (elliptical) galaxies, both BCGs and non-BCG ellipticals, that are three orders of magnitude more common than star-forming (spiral) galaxies (Davidzon et al., 2017), particularly in dense environments (e.g., Scoville et al., 2013; Darvish et al., 2016). By starting with this low-redshift population, and then identifying the progenitor populations responsible, we examine the role of internal versus external processes at a time when proto-clusters have not finished merging into the virialized clusters observed at present day. In situ star formation (SF) refers to any SF processes within a galaxy, such as those triggered through secular (disk instability; Livermore et al., 2012, 2015) or external stimuli (interactions or gas-rich mergers; Sanders & Mirabel, 1996; Kartaltepe et al., 2012). Ex situ stellar mass growth is the accretion of stellar mass formed outside of the galaxy, such as the old stellar populations composing the secondary galaxy in a galaxy merger (e.g., Pillepich et al., 2015).

The current hierarchical model of galaxy formation relies on a combination of stellar mass growth through internal processes and mergers that build galaxies from several individual galaxies. Mergers directly inject stellar mass into the primary galaxy concurrent with in situ SF triggered if either galaxy is gas rich. The relative roles of internal processes and mergers have been simulated at length in large-scale cosmological models (e.g., the Millennium-II Simulation (MS-II) & the Illustris project; Springel et al., 2005; Vogelsberger et al., 2014) and zoom-in models (Narayanan et al., 2010, 2015; Ragone-Figueroa et al., 2018). While there is still tension with observations of properties such as morphology (e.g., Genel et al., 2014), many predictions from simulations, such as merger rates, remain within observational errors (Lotz et al., 2011; Rodriguez-Gomez et al., 2017; Duncan et al., 2019). The total stellar mass growth of a BCG, or any massive galaxy, can be summarized by two drivers with multiple contributors.

In situ star formation

1. (a)

Self-sourced—a galaxy’s cold gas supply condenses via local or gravitational interaction/flyby. 2. (b)

ICM precipitation—cooling gas from the cluster ICM precipitates down the potential well and forms stars. 3. (c)

Wet mergers—stars form due to the shocks and greater gas supply provided during a gas-rich merger. 2. 2.

External delivery of stellar mass

1. (a)

Wet and dry mergers—immediate contribution of an existing, older stellar population. 2. (b)

Local interaction—capture of stellar mass from an interaction with an unbound secondary galaxy.

The relative role of each of these contributors is a rich field of research, with many works specializing in one or more facet of this picture (e.g., O’Dea et al., 2008; Naab et al., 2009; Loubser & Sánchez-Blázquez, 2012; Webb et al., 2015a, 2017; Contini et al., 2018). Another factor is the role of local environment on the stellar mass assembly processes described above. A correlation between local environment and SF is observed at low redshift (e.g., Oemler, 1976; Dressler, 1980; Kauffmann et al., 2004; Darvish et al., 2018), but is still a subject of debate at high redshift (e.g., Elbaz et al., 2007; Scoville et al., 2013; Darvish et al., 2016). Active galactic nucleus (AGN) activity as a function of environment at high redshift is likewise under investigation (e.g., Ehlert et al., 2013; Martini et al., 2013; Umehata et al., 2015; Alberts et al., 2016; Yang et al., 2018; Macuga et al., 2019).

In this paper, we seek to understand the roles of in situ SF and external delivery toward the total growth of BCG progenitors as a function of environment and redshift. However, identifying progenitors to specific populations in the present-day universe remains one of the most difficult, but critical, elements of galaxy evolution studies. Due to the long timescales involved, we require a method that connects galaxies across redshift and is independent of SFR characteristics such as color. This is often implemented as a stellar mass selection as a function of redshift and descendent stellar mass.

A relatively straightforward method for identifying the progenitors of massive galaxies is to assume that the most massive galaxies in clusters and groups remain the most massive throughout their history (van Dokkum et al., 2010). One implementation of this concept is the constant comoving number density method (Leja et al., 2013; Patel et al., 2013; van Dokkum et al., 2013), which assumes that the number of galaxies within a given comoving volume of the universe is constant. As the galaxies age within the volume, the population within each stellar mass bin will grow in stellar mass with time, and inhabit the new stellar mass bin at an equal density as their old bin. By selecting galaxies at high redshift that inhabit the universe at an equal cumulative number density, one selects the direct progenitors of the low-redshift population. This method can be corrected to include merger partners by increasing the cumulative comoving number density cutoff (i.e., lowering the stellar mass cutoff), and this evolving cutoff has been investigated using semi-empirical approaches or by using cosmological simulations (e.g., Marchesini et al., 2014; Torrey et al., 2015).

We use the evolving comoving cumulative number density selection method to identify the progenitors of BCGs out to redshift $z\sim 3$. We select a descendant stellar mass cutoff of 10${}^{11.5}\;M_{\odot}$ at $z<0.35$ to ensure we track progenitors to BCGs in the low-redshift universe. To characterize each progenitor, we fit the spectral energy distributions (SEDs) of each target from the far-ultraviolet (FUV) to far-infrared (FIR) to ensure sensitivity to low and obscured SFRs. Additionally, we build upon previous studies (e.g., Hill et al., 2017; Torrey et al., 2017) by correlating our progenitor characteristics to the local galaxy environment using redshift-sensitive adaptive kernel density maps of the COSMOS field (Darvish et al., 2015, 2017).

In Section 2 we discuss the data we used to construct the SEDs of our sample. In Section 3 we describe how we implemented an evolving comoving number density selection function to acquire our BCG progenitor sample. The median SED fits of our sample are discussed in Section 4. We consider the physical implications of our results in Section 5 and we summarize our findings in Section 6. We use the $\Lambda$CDM standard cosmological parameters of $H_{0}$ = 70 Mpc*-1* km s*-1*, $\Omega_{M}$ = 0.3, and $\Omega_{\mathrm{vac}}$ = 0.7 and a Chabrier (2003) initial mass function.

2. COSMOS Multiwavelength Data

We first require a large contiguous area and volume to examine progenitors across the full range of cosmological environments. Additionally, we need well-sampled SEDs from the FUV to the FIR to constrain the stellar mass and SFRs of the massive galaxies within these environments. These requirements are uniquely satisfied by the COSMOS field (Scoville et al., 2007). COSMOS is a 2 sq. deg. field centered at R.A.(J2000) = 10:00:28.600, decl.(J2000) = +02:12:21.00 observed from the X-ray to the radio. We use the public COSMOS photometric data catalog, COSMOS2015 (Laigle et al., 2016), to select progenitors and construct their SEDs from the FUV to FIR. Below, we concisely summarize each of the COSMOS datasets used in this study; for further details of COSMOS2015 see Laigle et al. (2016).

2.1. Photometric and Spectroscopic Redshifts

The quality of our SED fits and environment estimates depends on accurate distance measurements. To that end, we start with COSMOS photometric redshift measurements (Ilbert et al., 2008) in the COSMOS2015 catalog (Laigle et al., 2016). We then crossmatch the COSMOS2015 progenitors’ coordinates with the spectroscopic detections in the COSMOS spectroscopic catalog to the closest match within 1″  of the original COSMOS2015 location. Our BCG progenitor sample includes spectroscopic redshifts available from the 3D-Hubble Space Telescope (HST) Survey (Brammer et al., 2012; Momcheva et al., 2016), Keck DEIMOS 10K Spectroscopic Survey (Hasinger et al., 2018), FMOS-COSMOS Survey (Kashino et al., 2013; Kartaltepe et al., 2015b; Silverman et al., 2015), the Gemini GMOS-S spectra of Balogh et al. (2011), COSMOS AGN Spectroscopic Survey (Trump et al., 2007, 2009), the Keck LRIS spectra of Casey et al. (2017), MOIRCS Deep Survey (Yoshikawa et al., 2010), the Keck MOSFIRE spectra of Trakhtenbrot et al. (2016), the Keck DEIMOS spectra of Capak et al. (2011) and Mobasher (2016), the SINFONI spectra of Perna et al. (2015), the zCOSMOS Survey (Lilly et al., 2007), the VIMOS Ultra-Deep Survey (Le Fèvre et al., 2015), and the HST/WFC3 grim spectra of Krogager et al. (2014).

2.2. GALEX

We use FUV and near-ultraviolet (NUV) band point-spreadfunction (PSF)-fit photometric magnitudes from the Galaxy Evolution Explorer (GALEX; Martin et al., 2005) to constrain the degree of unobscured SF. GALEX observations in the COSMOS2015 catalog were originally reduced by Zamojski et al. (2007). We correct the GALEX FUV-Spitzer IRAC4 observations for Milky Way foreground extinction using a galactic reddening of $R_{v}$ = 3.1 (Morrissey et al., 2007) and E(B V) values from the dust maps of Schlegel et al. (1998).

2.3. Canada–France–Hawaii Telescope

Canada-France-Hawaii Telescope (CFHT)/MegaPrime (Aune et al., 2003; Boulade et al., 2003) $u^{*}$ observations of the COSMOS field were taken with a consistent depth of m${}_{u^{*}}$ $\sim$ 26.4 (Capak et al., 2007).

2.4. Subaru

To constrain our stellar mass estimates, we need sensitive optical continuum measurements across as many bands as possible. Therefore, we include Subaru/Suprime-Cam optical observations using five broadband filters (B, V, R, $i+$, $z++$) and 11 medium-band filters (IA427, IA464, IA484, IA505, IA527, IA574, IA624, IA679, IA738, IA767, and IA827), observed to a $3\sigma$ depth of $m_{AB}\sim 25.2$ or deeper. Out of these filters, IA464 suffers the worst resolution with a PSF FWHM of 1.89” (Taniguchi et al., 2007, 2015).

2.5. Vista

We include Vista/VIRCAM (Sutherland et al., 2015) J-, H-, and K-band observations from the UltraVISTA-DR2 survey (McCracken et al., 2012). These have limiting magnitudes $m_{AB}$ of 24.7, 24.3, and 24.0 respectively.

2.6. Spitzer Observations

Further IR observations of the COSMOS field are available from the SPLASH (Steinhardt et al., 2014) and S-COSMOS surveys (Sanders et al., 2007) using the Spitzer Space Telescope (Spitzer; Werner et al., 2004). Data from Spitzer’s Infrared Array Camera (IRAC)’s 3.6, 4.5, 5.7, and 8.0 µm  channels (for more information, see Fazio et al., 2004) have PSF widths 1$\farcs$6, 1$\farcs$6, 1$\farcs$8, and 1$\farcs$9, respectively, and were observed down to a 3$\sigma$ depth of $m_{AB}$ of 25.5, 25.5, 23.0, and 22.9, respectively. We also use MIPS (Rieke et al., 2004) 24 $\micron$ measurements, which act as an important constraint to the total IR fit. MIPS 24 $\micron$ was observed to a 5$\sigma$ depth of 71 $\mu$Jy.

2.7. Herschel Observations

To accurately estimate the obscured AGN and SF activity in our sample, we use FIR observations to constrain the FIR peak. COSMOS FIR observations include the Herschel Space Observatory (Herschel; Pilbratt et al., 2010) Photoconductor Array Camera and Spectrometer (PACS; Poglitsch et al., 2010) 100 and 160 µm  bands and Spectral and Photometric Imaging Receiver (SPIRE) 250, 350, and 500 µm  bands. The PACS Evolutionary Probe (Lutz et al., 2011) data used in this study were observed down to a 3$\sigma$ depth of 5 and 10.2 mJy for 100 and 160 µm  bands, respectively. Herschel Multi-tiered Extragalactic Survey (Oliver et al., 2012) SPIRE 250 µm, 350 µm, and 500µm  observations reach a 3$\sigma$ depth of 8.1, 10.7 mJy, 15.4 mJy, respectively.

Blending is possible at the longest wavelengths due to the large beam size (18.1″  at 250 µm, 24.9″  at 350 µm, 36.6″  at 500 µm), therefore de-blending was performed band by band using the next most-resolved observation starting with Spitzer IRAC observations as priors (Lee et al., 2010). Blending is the largest issue at 500 µm, so to ensure that the 500 $\micron$  detections correspond to our progenitor sample, we use the ‘clean index’ from Elbaz et al. (2011). In brief, a progenitor’s 500 $\micron$ data are rejected if there is more than one bright neighbor within 1.1x the FWHM at 24 $\micron$ or there are bright neighbors within 1.1$\times$ the FWHM of the respective Herschel FIR band, with a bright neighbor defined as $S_{\mathrm{Neighbor}}$/$S_{\mathrm{Target}}>$ 0.5.

3. Methods

3.1. Cumulative Number Density Selection Function

Due to the long timescales involved in an individual galaxy’s evolution, large samples over many epochs are needed to probe the growth of stellar mass in the universe. Selecting a population of progenitors at high redshift is difficult because the entire point of the exercise is that an important characteristic of the progenitor population remains unknown. The constant comoving cumulative number density method introduced by van Dokkum et al. (2010) was a solution to this problem. This method assumes a constant number of galaxies in the past and present universe. To identify the progenitors of massive galaxies in the present day, one would measure the density they inhabit in a volume element of the present-day universe and identify a corresponding population at high redshift that occupied the universe at the same volume density (illustrated in Figure 1).

This method has been successfully used to identify some of the progenitors for today’s population of massive galaxies (van Dokkum et al., 2010; Patel et al., 2013; Morishita et al., 2015), but it is biased toward higher-mass progenitors. This bias arises because the galaxies that merge onto the selected progenitor galaxy are not included. Secondary galaxies are often lower mass (of higher number density) and therefore their exclusion results in a bias toward lower number density and higher-mass progenitors (Mundy et al., 2015). Progenitors that evolved via a merger-dominant evolutionary path would not be selected, as these low-mass progenitors would inhabit a lower stellar mass (and higher number density) parameter space than the rest of the progenitor population. The disparity between predicted and ‘actual’ progenitor stellar masses could reach a factor of two at high redshift (Leja et al., 2013; Mundy et al., 2015). To better correct for the effect of mergers, an evolving number density method was developed by tracking the median progenitor mass in cosmological simulations (Mundy et al., 2015; Torrey et al., 2015; Jaacks et al., 2016).

An evolving number density selection is derived by identifying a population in a low-redshift slice of the universe and tracing back in redshift the number density that corresponds to the median stellar mass of their most massive progenitors. This method has been used with multiple simulations, and our analysis includes the Illustris simulation number density evolution tracks (Torrey et al., 2017, 2015), with comparison tracks from MS-II (Boylan-Kolchin et al., 2009) and a constant number density cut as discussed in Section 5.

The Illustris simulation is a cosmological hydrodynamic simulation of volume 1.2079 $\times$106 Mpc3 that includes dark and baryonic physics (Vogelsberger et al., 2014). Baryonic processes include radiative cooling, SF and supernova feedback, and AGN feedback. Baryonic mass is resolved with particles of mass 1.3 $\times$ 106 $M_{\odot}$. The inclusion of baryonic physics enables Illustris to more accurately model the evolution of SF at late times in the universe as AGN and outflows heat and push out the gas supply within each galaxy.

To estimate stellar masses of progenitors observed in COSMOS, we use Davidzon et al. (2017)’s SMF to identify the cumulative number density of galaxies with stellar masses above 1011.5 $M_{\odot}$ below $z\sim 0.35$, which corresponds to the regime of low-redshift BCGs (e.g., Lidman et al., 2012; Liu et al., 2012; Fraser-McKelvie et al., 2014; Cooke et al., 2016). This high stellar mass selection specifically examines high-mass BCGs commonly hosted in clusters with total mass $>10^{14}M_{\odot}$ (e.g., Cavagnolo et al., 2009; Postman et al., 2012; Lin et al., 2013) to better understand how the most massive galaxies acquired their stellar populations. We choose to use the SMF of the total population and not only the passive population, as the latter outnumbers the active population by a factor of ten at this mass and redshift and we do not wish to bias our selection toward the most passive progenitors. The results shown in Section 5 are consistent within errors if we exclude or include active descendants at low redshift. Once the density of the total descendant population is identified, we determine the cumulative comoving number density of BCG progenitors in redshift slices out to $z\sim 3$ using an Illustris-derived number density evolution relation (Equation (A4) of Torrey et al., 2017),

[TABLE]

[TABLE]

where $N_{0}$ represents the cumulative comoving number density of BCGs below $z\sim 0.35$, and the constants $A^{\prime}_{0,1,2}$ and $B^{\prime}_{0,1,2}$ are best-fit parameters to the backward median number density evolution.

3.2. Application to Observations

To construct our final sample, we take the cumulative comoving number density predictions estimated above and calculate the stellar mass corresponding to the progenitor cumulative number density in that redshift slice (red line in Figure 2) given the observed COSMOS SMF (Davidzon et al., 2017). To maintain sample completeness, we fit the SEDs of an expanded sample down to a lower stellar mass limit than required by our Illustris-derived selection function. By fitting COSMOS galaxies within a factor of 3 of the Illustris selection function, we hope to recover any galaxies with underestimated masses in their initial fit in the COSMOS2015 catalog. After fitting our initial sample of prospective progenitors (dotted blue line in Figure 2), we select galaxies that have a fit stellar mass (described in Section 4) above the Illustris evolving comoving number density method mass cut for that given redshift bin (red line in Figure 2). Our final sample includes 1444 BCG progenitors from $0.5<z<3.0$.

3.3. Comparison Sample Selection Methods

We also compare our results based on the Illustris selection function to other commonly used selection methods. The MS-II simulation (Boylan-Kolchin et al., 2009) is a dark-matter-only simulation and does not include baryonic physics. Number density is thus derived from dark matter halos characteristic of our low-redshift massive elliptical sample. For MS-II halos, we identify the stellar mass of a given halo mass’ s baryonic counterpart using Guo et al. (2013) Millennium WMAP 7 (MR7) semi-analytic galaxy models scaled to the MS-II simulation. We start similarly to our Illustris method by using Davidzon et al. (2017)’s SMF to identify the cumulative number density of galaxies with stellar masses above 1011.5 M*⊙* at $z\sim 0.35$. Then we identify the dark matter halo mass corresponding to this number density cutoff. We identify the most massive halo progenitor in each higher-redshift bin, and find the number density corresponding to the median of the most massive halo. This number density is then used with the corresponding Davidzon et al. (2017) SMF at that redshift to identify the galaxy progenitor population from Guo et al. (2013)’s Millennium-II mock galaxy catalog. We make the assumption that a single BCG progenitor lies in each halo. The MS-II derived evolution track is shown as the teal solid line in Figure 2. We discuss how the MS-II-selected sample results compare to our primary Illustris-selected results in Section 5.6.

Finally, we wish to compare our results to those derived from the original constant number density method. Our constant number density comparison sample is chosen by identifying the cumulative comoving number density of ellipticals with stellar masses 1011.5 $M_{\odot}$ at $z\sim 0.35$ in Davidzon et al. (2017)’s SMF, and identifying the stellar mass of the cumulative population above that density in all redshift bins out to $z\sim 3$. This stellar mass and redshift combination corresponds to the massive BCG population in the low- to intermediate-redshift universe (e.g., Lidman et al., 2012; Liu et al., 2012; Fraser-McKelvie et al., 2014; Cooke et al., 2016). Our constant number density evolution track is shown as the orange solid line in Figure 2. A comparison between the constant number density selection method and the evolving number density method used for the main results of this work is included in Section 5.7.

3.4. Quantifying Local Environment

With progenitor stellar mass selection criteria calculated, our next step is to select galaxies likely to merge onto the BCG ‘seeds’ over time given their local environment. The sample selection discussed in Section 3.1 identified progenitors, irrespective of environment, that are predicted to relocate to the overdense regions where present day BCGs ultimately reside. Therefore we use density maps of the COSMOS field to identify progenitors already in overdense regions (very likely to merge into a BCG within a Hubble time) and underdense regions (less likely to merge into a BCG).

We quantify local environment using density field maps of the COSMOS field produced using the ‘weighted’ adaptive kernel smoothing approach of Darvish et al. (2015). The density field is evaluated in a series of overlapping redshift slices. Darvish et al. (2015) assign galaxy weights that are proportional to the likelihood of a galaxy belonging to a redshift. They are estimated by measuring what fraction of a galaxy’s redshift distribution is within a particular redshift slice. Given these weights, the density field is then adaptively smoothed with a Gaussian kernel whose ‘global’ width is 0.5 Mpc. The adaptive widths are evaluated around this global width, depending on how sparse or dense the local neighborhood of each galaxy is (see Darvish et al., 2015, for details).

3.5. Identification of AGN

We also quantify the contribution of the AGN emission to the IR regime of progenitor SEDs. This step is important, as a dominant AGN could cause the SED fit solution to overestimate stellar mass, and an IR-bright AGN could result in a fit with an overestimated SFR. We crossmatch our sample with the Chandra X-ray observations of the COSMOS field (Elvis et al., 2009; Civano et al., 2016; Marchesi et al., 2016; Lanzuisi et al., 2017) to identify X-ray-bright AGN. We also use the Spitzer IRAC photometry available for the COSMOS field to identify obscured AGN through their mid-infrared (MIR) colors using the Donley et al. (2012) IRAC color criteria. Radio AGN are identified by crossmatching our progenitor sample with the radio AGN identified in the VLA-COSMOS 3 GHz Large Project (Delvecchio et al., 2017; Smolčić et al., 2017, 2008). Our full progenitor and AGN sample selection totals and statistics are listed in Table 1.

4. SED Fitting

4.1. Software and Methodology

To determine the SFR and total stellar mass of each progenitor, we require a way to reliably fit thousands of galaxy SEDs from rest-frame FUV to FIR. To do so, we choose the SED fitting software MAGPHYS (da Cunha et al., 2008) with 50,000 stellar population models generated from Bruzual & Charlot (2003) and 50,000 infrared dust models from da Cunha et al. (2008).

MAGPHYS fits the SED of each galaxy in two stages. First, it estimates the best-fit stellar model between $0.0912\micron\lesssim\lambda\lesssim 10\micron$ from the Bruzual & Charlot (2003) stellar population models while taking into account dust obscuration. Once the optical component is fit, MAGPHYS fits the infrared from 2.5 to 1000 µm  while taking the total luminosity of obscured starlight found during the optical fitting as a prior. FIR characteristics are calculated between 8 and 1000 µm. This ensures energy balance between the energy absorbed and emitted by the dust component.

At time of writing, the public version of MAGPHYS does not include an AGN component. This is an important consideration, as additional luminosity across the optical and near-IR from an AGN accretion disk and the dusty torus surrounding it may skew stellar mass estimates to higher values. Additionally, any UV or IR emission may also cause SFRs to be overestimated. Therefore to identify potential overestimations due to a lack of an AGN fitting component, we investigated several goodness-of-fit metrics to identify which progenitors require an AGN component. We find that a fractional difference cutoff of 0.4 between the MAGPHYS fit 8.0 $\micron$ (Spitzer IRAC4) flux and the observed value worked most consistently to identify sources that require an AGN component in the MIR (for further details, see K.D. Tyler et al. 2019, in preparation). This parameter measures how poorly the stellar and dust models fit the slope of the MIR, which could indicate the requirement of a third (AGN) component. We find that $\sim$1% of our progenitor sample fits have MAGPHYS fractional 8.0 $\micron$ residuals $>0.4$, and we refit this subsample using the package SED3FIT. By choosing to select and fit AGN components to only the galaxies that have poor IR fits without an AGN, we prevent the overfitting of our sample and the possible underestimation of SFRs that an overestimation of AGN activity and frequency would induce.

SED3FIT (Berta et al., 2013) works similarly to MAGPHYS by fitting a series of optical models first, and then using the optical results as a prior during the fitting of the IR component. SED3FIT, however, includes AGN torus and accretion disk emission (Feltre et al., 2012) across both components while maintaining energy balance. For the IR-bright AGN identified through their $8\mu$m residuals, we use the SFR and stellar mass estimates from SED3FIT instead of MAGPHYS in the following results.

4.2. AGN Template Library

The original AGN template library included in the public SED3FIT distribution included 10 AGN templates from the Feltre et al. (2012) AGN template library. These templates were selected at extreme viewing angles, effectively probing archetypical Type-1 and Type-2 AGN. These models assume either negligible or strong extinction with respect to line of sight. While satisfactory for first-order fits, the model of the obscuring torus surrounding the supermassive black hole in an AGN has changed in recent years to a clumpy torus with regions of high and low obscuration at high and low viewing angles (e.g., Krolik & Begelman, 1988; Shi et al., 2006; Markowitz et al., 2014). Therefore, we expand our AGN template library to include 240 models (K.D. Tyler et al. 2019, in preparation) equally probing between 0$\degree$ and 90$\degree$. Shown in Figure 3, this new library contains models that fit AGN at intermediate obscurations and viewing angles.

In Figure 4, we plot the median SED of the progenitors found in three environments in six redshift bins with $1\sigma$ confidence intervals (shaded regions). We classify each environment as dense ($\frac{\delta}{\delta_{median}}>$ 2), intergroup (1 $<$ $\frac{\delta}{\delta_{median}}<$ 2), or field ($\frac{\delta}{\delta_{median}}<$ 1). In the top panel, we find that progenitors in dense regions have SEDs similar to those in other environments at high redshift ($z>1.1$), but have $\sim\frac{1}{5}$ the MIR emission (10–20$\mu$m) at $z<1.1$ in comparison to progenitors at lower density at low redshift. This indicates a divergence in evolution due to environment at low redshift, where massive progenitors in high-density environments undergo a more efficient removal of their warm gas and dust supplies.

5. Results and Discussion

To characterize our total sample of 1444 BCG progenitors, we use the stellar parameters fit by MAGPHYS for 1407 progenitors and the estimates from SED3FIT for 37 progenitors due to their significant 8.0 $\micron$ residuals in their original MAGPHYS fits (hereafter referred to as SED3FIT AGN). The SED3FIT AGN fraction is $<1\%$ in any redshift bin, with an overall AGN fraction found through X-ray, MIR, and radio methods of $\sim 8\%$ (Table 1). For comparison, previous Sunyaev-Zeldovich-selected studies have found X-ray cavities indicating past AGN activity in $\sim$7$\%$ of low-redshift BCGs (Hlavacek-Larrondo et al., 2015).

5.1. SF and Stellar Mass Estimates

Our sample spans an order of magnitude in stellar mass from 1011 – 1012 $M_{\odot}$ at $z\sim 3$ (Figure 5), motivated by the stellar mass selection function shown in Figure 2. Our sample below $z\sim 0.8$ is subject to small-number statistics due to the low volume available in COSMOS at this distance ($\sim$ 4.4 $\times$ 106 Mpc3).

Shown in Figure 5, both SFR (averaged over 108 yr) and specific SFR (sSFR) span three orders of magnitude from $0.5<z<3$. Progenitors do not span this parameter space equally, with a concentration of highly star forming progenitors with sSFR $\sim$ 10*-9.25* yr*-1* at $z>2$, that transition to a lower sSFR of 10*-10.75* yr*-1* at $z<2$. This is partly driven by our stellar mass selection function, which is limiting the stellar mass parameter space to higher stellar masses at low redshift. Our observed relation at the bottom of Figure 5 is systematically an order of magnitude below the sSFR–$z$ correlation found for less-massive (1010.5 $M_{\odot}$) galaxies found in the COSMOS field (Davidzon et al., 2018), consistent with galactic ‘downsizing’ where less-massive galaxies continue to form stars while massive galaxies are quenching (Cowie et al., 1996; Pérez‐-González et al., 2008). This will yield descendants dominated by the old stellar populations, just as observed in BCGs in the nearby universe (Loubser et al., 2009).

5.2. Evolution of Stellar Mass Growth

Our goal is to investigate what role in situ SF plays in the total growth of stellar mass in BCG progenitors. To do so, we need to compare the growth solely due to SF to a measurement of total stellar mass growth. We define our total growth rate to be the change in median stellar mass between each redshift bin in units of $M_{\odot}$ yr*-1*. We then perform a polynomial fit to these growth rates at each bin boundary and plot the result as the solid green line in Figure 6.

We determine the growth rate of each progenitor due to in situ SF by taking the SFRs from our SED fits and applying a reduction of 50$\%$. This scale factor used by Brinchmann et al. (2004) and van Dokkum et al. (2010) corrects for stellar mass loss over a Gyr timescale. This correction to SFR produces a stellar mass growth measurement that considers the net growth influenced by SF and star destruction processes. Stellar mass loss occurs as the population of stars born during an SF episode lose their massive members over time as they age off the main sequence. The final in situ growth rate represents the total stellar mass per year generated in a galaxy that will remain on multi-Gyr timescales. After we apply this correction, we calculate the in situ growth rate for progenitors in dense ( $\frac{\delta}{\delta_{\mathrm{median}}}>$ 2 ) and field environments ($\frac{\delta}{\delta_{median}}<$ 1 ) environments and plot the mean stellar mass growth rates due to star formation for each sample as red and blue lines, respectively, in Figure 6. $\delta$ is defined as the density at that location while $\delta_{\mathrm{median}}$ is the median density of the redshift slice.

Finally, we wish to understand how the stellar mass directly delivered via gas-rich and gas-poor mergers compares to the above values. Therefore we integrate the stellar mass deposition rate ($M_{\odot}$ yr*-1*) from mergers of secondary galaxies with mass ratios of 1:1 to 1:10*-5* of the progenitor’s stellar mass in that redshift bin using the merger rates provided by the Illustris simulation (Rodriguez-Gomez et al., 2017). This is plotted as the dashed gray line in Figure 6.

Following Hill et al. (2017), we plot the above values as $dM/dt$ ($M_{\odot}$ yr*-1*) versus redshift from SF, mergers, and all sources in Figure 6. From this comparison, we identify three epochs of BCG progenitor growth.

5.2.1 The SF-Dominated Epoch ($z>\mathit{2.25}$)

At $z>2.25$, the in situ SF estimated in progenitors in all environments is consistent with the total growth rate of the median stellar mass of our sample. This consistency indicates that progenitor growth is dominated by active SF with contributions of direct stellar mass injection and gas via gas-rich and gas-poor mergers. The total stellar mass growth rate shown in Figure 6 is motivated by the stellar mass selection function that was used to originally select our sample; however, we may still conclude that the active SFRs via secular or wet merging processes are high enough to be the primary method of individual progenitor growth. Dry mergers are unnecessary, but may still contribute to the highest-mass cases.

5.2.2 The Transitionary Epoch ($\mathit{1.25}<z<\mathit{2.25}$)

Between $1.25<z<2.25$, progenitors are forming their in situ stellar mass at a rate $\sim$60–75$\%$ of the observed median total mass growth rate. In this era, progenitors are still forming their mass predominantly through in situ SF independent of environment; however, additional methods of stellar mass generation are required. It is during this time that the total stellar mass delivered by mergers becomes important to BCG progenitor evolution, as the mass formed via SF (red/blue lines in Figure 6) is insufficient to account for the observed total stellar mass growth rate (green line in Figure 6).

5.2.3 The Dry Merger Epoch ($z<\mathit{1.25}$)

Finally, at $z<1.25$ the in situ SF of BCG progenitors in dense environments sharply declines with redshift, while massive field galaxies maintain SF activity down to lower-redshifts. In our lowest redshift bin ($z=0.55$), the limited area of the COSMOS field introduces small-number statistics. In this bin, we observe only one galaxy above our mass cut in a field environment, compared to 16 in dense environments. We see a general trend where SF is more than an order of magnitude insufficient to account for the total stellar mass growth rate. This era requires ex situ mass delivery systems, specifically dry mergers, to be the dominant growth mechanism.

5.3. The Onset of Quenching

We contextualize the diminishment of in situ SF toward progenitor growth by examining how rest-frame colors and morphological parameters such as the Sérsic index evolve with redshift. Shown in Figure 7, we plot the rest-frame NUV-r and r–J colors (Ilbert et al., 2013) of our progenitor sample in each of our six redshift bins with Sérsic index from the COSMOS Zurich Morphological Catalog (Sargent et al., 2007; Scarlata et al., 2007) as the colorbar. We use the rest-frame absolute magnitudes provided by the COSMOS2015 catalog. Rest-frame magnitudes were calculated using the observed nearest filter magnitude and a $k$-correction estimated from the SED fitting performed by Laigle et al. (2016). Our SFR estimates are consistent with their target’s NUVrJ classification, with the star-forming galaxy subpopulation hosting median SFRs 10 times the SFR of the quiescent population in our lowest-redshift bin, and five times the SFR of the quiescent population in our highest-redshift bin. A caveat remains that this diagnostic may misclassify star-forming galaxies as quiescent, and contaminate $>10\%$ of quiescent selected galaxies (Ilbert et al., 2013). Therefore we compare the total number of galaxies classified as quiescent to the number of quiescent galaxies above the star-forming main sequence, the correlation between SF and stellar mass indicating the median star-forming activity for a galaxy of a given stellar mass and redshift (e.g., Schreiber et al., 2015). We find that $\sim$3-5% of progenitors classified as quiescent lie above the star-forming main sequence, potentially the misclassified star-forming galaxies found in Ilbert et al. (2013). This contamination rate is low enough that the results shown in Figure 7 are not affected.

Sérsic indices were estimated using the HST/ACS F814W observations of the COSMOS field. NUVrJ is more sensitive to recent SF evolution than UVJ-derived color measurements (Martin et al., 2007; Davidzon et al., 2017), an important consideration in determining when a galaxy has begun quenching. Position in this color space has also been shown to correlate with sSFR as the galaxy’s star forming efficiency changes (Martin et al., 2007; Arnouts et al., 2013; Ilbert et al., 2015). We find that the stellar populations of our BCG progenitors transition to a quiescent state between our $0.8<z<1.1$ and $1.1<z<1.5$ redshift bins, including those without morphological classifications. This is consistent with our results from Figure 6, where in situ SF becomes a negligible contributor to a BCG’s evolution between $1.0<z<1.3$. Cotemporal with the transition to a quiescent stellar population, the morphology of the progenitor population also changes. The median Sérsic index in each redshift bin (Figure 7) is constant within the scatter, centered at $\langle$n$\rangle\sim 2.5$. These median Sérsic indices are consistent with those found for S0 galaxies in cluster environments (e.g., D’Onofrio et al., 2015), suggesting that galaxies in our progenitor sample have hosted a composite bulge+disk structure since $z\sim 2$.

5.4. Visual Classification of Morphology

To provide context to the light profile evolution observed in Section 5.3, we classify the morphological state of all progenitors imaged in the COSMOS field ACS mosaic (Koekemoer et al., 2007), totaling 1372 of our original sample of 1444. Progenitors not included are outside the bounds of the ACS mosaic. The COSMOS ACS mosaic includes HST/ACS F814W broadband filter observations with a final pixel scale of 0.03″ per pixel and a 3$\sigma$ surface brightness density of $m_{AB}$ $\sim$ 25.1 mag arcsec*-2* (Wen & Zheng, 2016). This spatial resolution corresponds to 183 pc per pixel at $z\sim 0.5$ and 231 pc per pixel at $z\sim 3$. We also supplement this observation set with any publicly available images from HST/WFC3 F105W, F110W, F125W, F140W, and F160W. This subset of 381 targets is most useful in our high-redshift bins as the peak of the stellar continuum is redshifted into the WFC3 IR broad bands. The following results shown in this work are consistent independent of the inclusion or exclusion of targets with WFC3 observations. To ensure consistency of classification, we utilize the morphological classifications scheme and user interface implemented by the CANDELS team (Kartaltepe et al., 2015a; Kocevski, 2015) to classify each galaxy’s morphology (spheroid, disk, irregular) as well as interaction class (merger, interaction, non-interacting companion).

The user examines the HST image and labels the progenitor with as many flags as necessary to classify it. Our sample was classified by three authors (K.C.C., J.S.K., K.D.T.), and two of three authors must agree on a designation to assign a galaxy to any of the classification bins. For galaxies that lack the organized, symmetric morphology of a disk and/or spheroid, the irregular flag is used. To provide context to the primary classification, interaction classes are also used for the apparent stage of an interaction. The merger flag is used in highly disturbed cases where the primary and secondary galaxies of the mergers are nearly indistinguishable. The interaction flag is used for disturbances such as tidal tails where the interacting galaxies remain independently resolved. For more details on this classification scheme see Kartaltepe et al. (2015a). The distribution of morphologies for our BCG progenitor sample is shown in Figure 8. Our sample is predominantly spheroidal at low redshift, with a growing fraction of unclassifiable progenitors toward high redshift, where the rest-frame optical stellar continuum is redshifted out of the observed band. We see the majority of progenitors are spheroidal, indicating that a significant fraction of our total progenitor sample has already developed the spheroidal structure exhibited by present day BCGs. This appears to be consistent with recent work that identifies a population of massive spheroidals at high redshift (e.g., Cassata et al., 2013).

To compare the environmental effects discussed in Sections 5.2.1–5.2.3, we plot the net classification, or ‘diskiness’ of progenitors in our three environment bins (top half of Figure 9). Our diskiness parameter is defined by the averaging of disk flags (valued as +1) and spheroid flags (valued as -1). For example, a target with two spheroid flags and one disk flag will have a diskiness value of -0.33. Flags indicating disk- or bulge-dominated features are also included as values of $\pm$0.5, respectively. These flags are used to refine the classification of a galaxy with multiple primary flags and identify objects such as bulge-dominated S0 galaxies. We find that low-redshift progenitors have a decreasing likelihood of being labeled a disk with increasing local density. This trend is muddled in our highest-redshift bins, where we do not observe evidence for an environmental dependence. Next, we plot the diskiness parameters of progenitors identified as star-forming and quiescent through their NUVrJ colors (see Section 5.3) in the bottom half of Figure 9. We find a clear offset of average net classification for these two subsamples, where quiescent galaxies are more spheroidal on average than their star-forming counterparts below $z\sim 2.25$.

We next consider how often progenitors are classified as an irregular with redshift, and how this flag is used in combination with other flags such as disk or spheroid. We include example HST/ACS F814W images of the possible irregular combinations in Figure 10. In Figure 11, we plot the fraction of the total sample that has at least two irregular flags per target, identical to the confident irregular subsample in Figure 9. This total irregular sample is divided into six redshift bins with subsamples that include at least one spheroid, disk, or composite (spheroid+disk) flag. As these are confident irregulars (two or more irregular flags), we lower the threshold of additional flags to one to examine how often other characteristics such as disks or spheroids are assigned to this population. The frequency of irregulars in our sample declines with decreasing redshift. Also with redshift, we see the irregular population transition from a diverse range of morphologies to an irregular-spheroid and irregular-composite population toward present day. One result of interest is the lowest-redshift bin, in which disk components are identified in all the irregular candidates, a rare feature ($\leq$ 5%) in the total sample at $z<1$.

Finally, we examine the distribution of interaction classifications (Figure 12). To classify the state of interactions, we use flags to identify any ongoing ‘mergers’ in which the structure of both galaxies have been disturbed and combined. An ‘interaction’ flag is used when an interacting pair, identified using tidal features, are visibly distinct from one another. We note that these classifications do not distinguish between major and minor merger pairs. We observe a slowly decreasing interaction fraction with redshift ($\sim 25-15\%$), with a roughly constant merger fraction of $5\%$ of the total sample. We find a $1\sigma$ Poisson confidence interval of $\pm 5\%$ between $0.5<z<2.0$ and $\pm 3\%$ between $2.0<z<3.0$. Our results are consistent within errors of those calculated by Man et al. (2016) for the general galaxy population; they find a total interaction frequency of $\sim$20% at $z\sim 1$ using a lower-mass sample of pairs down to $10^{10.8}$ $M_{\odot}$. The lack of observed mergers in our lowest-redshift bin is consistent with a low-redshift, minor merger-dominated model for BCG growth (e.g., Edwards & Patton, 2012). It is also consistent with the results from Duncan et al. (2019), who found a merger fraction of $<10\%$ for their massive subsample (1010.3 $M_{\odot}$) in the CANDELS fields. The lack of interactions in our highest-redshift bin may be due to the difficulty in visually identifying low surface brightness features present in interactions at $z>2$. Overall this slowly evolving fraction of interactions is consistent with a model with fewer mergers at later times.

5.5. Comparison to Literature

Several works on massive galaxy evolution and the comparison of progenitor selection methods have been published in recent years (e.g., Zhao et al., 2016; Hill et al., 2017; Torrey et al., 2017), with results consistent with ours. Hill et al. (2017) selected massive galaxy progenitors using the number density evolution functions derived from the Bolshoi simulation (Klypin et al., 2011). Due to the dark-matter-only nature of the Bolshoi simulation, a comparison between our work and that done by Hill et al. (2017) helps illustrate how stellar mass selection functions derived using different combinations of dark and baryonic physics affect measured evolution.

Hill et al. (2017) used a lower stellar mass progenitor selection than ours. Since lower-mass galaxies build their stellar mass over longer timescale on average (mass downsizing; Cowie et al., 1996; Pérez‐-González et al., 2008), Hill et al. (2017)’s selection led to an in situ SF-dominated phase which persisted to lower redshift than ours, until $z\sim 1.75$. However, they found a fast transition to a quenched state from $z<1.75$ (Figure 1 of Hill et al., 2017). Dark-matter-only simulations producing a lower-mass selection function is hypothesized to be due to the lack of baryonic physics feedback effects that limit the speed of stellar mass growth, enabling the inclusion of faster evolving, lower-mass progenitors.

The median stellar mass of our sample increases by a factor of $\sim 80\%$ over $0.35<z<2.00$, and a factor of two when extrapolated to the present day from $z=2$. This rate of BCG stellar mass growth is slower than the factors of two from $z=1$ shown in BCGs residing in clusters selected from Spitzer observations (Lidman et al., 2012; Lin et al., 2013; Shankar et al., 2015). However, our results are in rough agreement with Zhao et al. (2016), who found that BCG progenitors grow by a factor of 2 from $z\sim 2$ to present day. The Zhao et al. (2016) sample was selected using a hybrid selection function which selected the most massive galaxies in the 38 most overdense regions identified in the CANDELS UDS field. The Zhao et al. (2016) sample also exhibited dry merger-dominated growth at $z<1$. Moreover, our work is consistent with additional works that indicate a slow rate of growth at low redshift (e.g., Tonini et al., 2012; Hill et al., 2017), shown in BCGs selected via the Maraston (2005) stellar population models or comoving number density cuts. Our observation of a transition to a dry merger-dominated epoch below $z\sim 1$ is also consistent with the Spitzer Wide-Area Infrared Extragalactic (SWIRE) Survey (Lonsdale et al., 2003) samples of Webb et al. (2015b). However, above $z\sim 1$, our sample includes progenitors with SFRs an order of magnitude lower than those of Webb et al. (2015b). This is believed to be primarily due to the difference in depth of the 24$\mu$m imaging in COSMOS versus SWIRE (71 $\mu$Jy versus $\sim$150 $\mu$Jy, respectively).

5.6. Comparison to the Millennium-II Simulation

As there is not yet a perfect cosmological simulation of the universe, it is worthwhile to examine how choosing a different simulation may change our results. Dark-matter-only simulations such as MS-II and Bolshoi, for example, have a different merger tree history than Illustris, which changes the rate ex situ stellar mass is accumulated. To test this, we retrieve the median cumulative number density of progenitor halos from MS-II (Boylan-Kolchin et al., 2009). To identify BCG progenitors simulated in MS-II, we identify the halo population inhabiting a comoving volume at the same density as observed low-redshift BCGs (Davidzon et al., 2017). We then identify the most massive progenitor halo corresponding to each descendant halo in ascending redshift slices to $z\sim 3$ and measure the cumulative comoving number density corresponding to the median most massive halo in that redshift slice. This number density is used to select mock galaxies from Guo et al. (2013), which represent the baryonic content of the MS-II BCG progenitor halos. The median Guo et al. (2013) sample stellar mass is used as to select observed galaxies in COSMOS based on their stellar mass (teal line in Figure 2).

We fit the SEDs of MS-II-selected progenitors identified in the COSMOS2015 catalog in the same manner as the Illustris-selected progenitors, as described in Section 3. Shown in Figure 13, the in situ growth of massive galaxy progenitors due to star formation is at least half of the total growth until $z\sim 1.5$, lower than our findings using Illustris-selected progenitors. However, the total growth rate track has a much higher normalization than Illustris, indicating a later, faster growth history in MS-II. To be physically consistent between SFR and total growth in MS-II selected progenitors, a much richer merger history is required to deliver mass onto the simulated BCG seeds in addition to their intrinsic SF. These differences between Illustris and MS-II may have several causes. One is the difference in feedback physics between the MS-II (Guo et al., 2013) semi-analytical models assigned to their halos and the baryonic physics in Illustris. A less efficient feedback process in MS-II would prevent the slow- or shut-down of highly star forming progenitors. This faster evolution produces a mass selection function with a larger slope, including galaxies with correspondingly lower SFRs.

5.7. Comparison to the Constant Density Method

We also compare our results to the constant cumulative number density selection method (van Dokkum et al., 2010). Following the constant number density selection function plotted in Figure 2, only the most massive progenitors are selected. As discussed in Section 3, a constant number density selection only selects progenitors that inhabit the universe at an equal density, an assumption that the number of galaxies in the universe is constant. As this ignores the effect of mergers, the median progenitor mass is not pulled down by the lower-mass secondaries. We use the same fits described in Section 3 to estimate the average stellar mass accretion via in situ SF in comparison with the evolution of the median stellar mass of the progenitor sample in Figure 14. We find that the combination of slow predicted stellar mass evolution and selection of very massive galaxies that are likely to be star forming result in a model where SF is dominant until $z<1$, after which gas poor mergers are necessary.

6. Summary and Conclusions

We examine the in situ SF and ex situ stellar deposition components of the total stellar mass growth of progenitors for the present day BCG population. By using FUV–FIR photometric data in the COSMOS field, we fit the SEDs for several samples of BCG progenitors and determine their SFR, sSFR, and stellar mass. Once corrected for stellar mass loss, we identify multiple epochs of BCG evolution based on the prominence, or lack thereof, of in situ SF toward stellar growth. Our estimates of in situ SF include episodes fueled by a progenitor’s own gas supply as well as those triggered by gas-rich mergers. Total stellar mass growth also includes the direct contribution of stellar mass by dry mergers. These behaviors are then compared to local environment to test for the presence of any dependence. Our results are as follows.

• •

BCG progenitor growth is dominated by in situ SF, either secular or merger driven, until approximately $z\sim 2.25$, after which galaxy mergers (gas rich and gas poor) contribute increasing fractions of ex situ stellar mass growth with time until $z\sim 1.25$. After this time, dry mergers become the dominant growth mechanism.

• •

Any dependence on local environment at $z>1$ is indistinguishable within our errors. The SFR–local density correlation observed today is observed in our $z\sim 0.8$ bin; however, volume effects limit our ability to measure this effect at lower redshifts. We also find that spheroidal progenitors become more common in dense environments ($\frac{\delta}{\delta_{median}}>$ 2) than field environments ($\frac{\delta}{\delta_{median}}<$ 1), which most often host spheroid+disk composites.

• •

The era of quenching, either through merger or AGN activity, is encountered by most progenitors by $z\sim 1.00-1.25$. This era is characterized by a shut-down of in situ SF and the transition from a disk+spheroid morphology to an overall spheroidal type.

• •

Across all redshifts, progenitors requiring an AGN-fitting component represent a small fraction ($\sim$ 1%) of the total progenitor population. SED3FIT AGN are most often observed at $z>1.5$, with 29/37 of our SED3FIT AGN sample observed at high redshift.

• •

Our progenitor sample consists of spheroid+disk composites until $z\sim 1.3$, after which we identify an environmental dependence on morphology. Progenitors in dense environments are labeled with a more spheroidal net classification as redshift decreases. This is primarily driven by the change in diskiness of the quiescent population in dense regions. Progenitors in intergroup regions ($1<\frac{\delta}{\delta_{median}}<$ 2) retain their composite structure. The few galaxies in field environments that meet our stellar mass cut increase in diskiness. Mergers and interactions are no more common than the general population ($<25\%$), with irregular galaxies exhibiting a greater range of morphological classifications with higher redshift.

We thank the referee for the constructive feedback that improved the clarity of this work. Support for this work was provided by NASA through grants HST-GO-13657.010-A and HST-AR-14298.004-A awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555. Support and observations made with the NASA/ESA Hubble Space Telescope, obtained from the Data Archive at the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., were also provided by NASA through grant NNX16AB36G as part of the Astrophysics Data Analysis Program. Spectral energy distribution fitting was performed using the scientific computing resources of the Rochester Institute of Technology ION cluster. B.D. acknowledges financial support from NASA through the Astrophysics Data Analysis Program (ADAP), grant No. NNX12AE20G, and the National Science Foundation, grant No. 1716907. We thank Paul Torrey for his council on evolving number density methodology. We thank the generous advice of Elisabete da Cunha and Stefano Berta on the use of MAGPHYS and SED3FIT, respectively. This research made use of Astropy,111http://www.astropy.org a community-developed core Python package for Astronomy (Robitaille et al., 2013; Price-Whelan et al., 2018). This research made use of the iPython environment (Perez & Granger, 2007) and the Python packages SciPy (Jones et al., 2001), NumPy (van der Walt et al., 2011), and Matplotlib (Hunter, 2007). Based on observations made with the NASA Galaxy Evolution Explorer. GALEX is operated for NASA by the California Institute of Technology under NASA contract NAS5-98034. MegaPrime/MegaCam, a joint project of CFHT and CEA/DAPNIA, at the Canada-France-Hawaii Telescope (CFHT) which is operated by the National Research Council (NRC) of Canada, the Institut National des Sciences de l’Univers of the Centre National de la Recherche Scientifique of France, and the University of Hawaii. Based [in part] on data collected at Subaru Telescope, which is operated by the National Astronomical Observatory of Japan. Based on observations obtained with WIRCam, a joint project of CFHT, Taiwan, Korea, Canada, France. This work is based [in part] on observations made with the Spitzer Space Telescope, which is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA. Herschel is an ESA space observatory with science instruments provided by European-led Principal Investigator consortia and with important participation from NASA. This work was supported by a NASA Keck PI Data Award, administered by the NASA Exoplanet Science Institute. Data presented herein were obtained at the W. M. Keck Observatory from telescope time allocated to the National Aeronautics and Space Administration through the agency’s scientific partnership with the California Institute of Technology and the University of California. The Observatory was made possible by the generous financial support of the W. M. Keck Foundation. The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Mauna Kea has always had within the indigenous Hawaiian community. We are most fortunate to have the opportunity to conduct observations from this mountain.