We present a robust method to constrain average galaxy star formation rates SFRs star formation histories SFHs and the intracluster light ICL as a function of halo mass. Our results are consistent with observed galaxy stellar mass functions specific star formation rates SSFRs and cosmic star formation rates CSFRs from z  0 to z  8. We consider the effects of a wide range of uncertainties on our results including those affecting stellar masses SFRs and the halo mass function at the heart of our analysis. As they are relevant to our method we also present new calibrations of the dark matter halo mass function halo mass accretion histories and halosubhalo merger rates out to z  8. We also provide new compilations of CSFRs and SSFRs more recent measurements are now consistent with the buildup of the cosmic stellar mass density at all redshifts. Implications of our work include halos near 10SUP12SUP M SUBSUB are the most efficient at forming stars at all redshifts the baryon conversion efficiency of massive halos drops markedly after z  2.5 consistent with theories of coldmode accretion the ICL for massive galaxies is expected to be significant out to at least z  11.5 and dwarf galaxies at low redshifts have higher stellar mass to halo mass ratios than previous expectations and form later than in most theoretical models. Finally we provide new fitting formulae for SFHs that are more accurate than the standard declining tau model. Our approach places a wide variety of observations relating to the SFH of galaxies into a selfconsistent framework based on the modern understanding of structure formation in CDM. Constraints on the stellar masshalo mass relationship and SFRs are available for download online.
['dark matter', 'galaxies: abundances', 'galaxies: evolution', 'methods: numerical', 'Astrophysics - Cosmology and Extragalactic Astrophysics']
The Average Star Formation Histories of Galaxies in Dark Matter Halos from z  08
{'THE AVERAGE STAR FORMATION HISTORIES OF GALAXIES IN DARK MATTER HALOS FROM Z = 0 -8': 'PETER S. BEHROOZI, RISA H. WECHSLER  \nKavli Institute for Particle Astrophysics and Cosmology; Physics Department, Stanford University; Department of Particle Physics and Astrophysics, SLAC National Accelerator Laboratory; Stanford, CA 94305', 'CHARLIE CONROY': 'Department of Astronomy & Astrophysics, University of California at Santa Cruz, Santa Cruz, CA 95064 Draft version November 27, 2024', 'ABSTRACT': 'We present a robust method to constrain average galaxy star formation rates, star formation histories, and the intracluster light as a function of halo mass. Our results are consistent with observed galaxy stellar mass functions, specific star formation rates, and cosmic star formation rates from z = 0 to z = 8. We consider the effects of a wide range of uncertainties on our results, including those affecting stellar masses, star formation rates, and the halo mass function at the heart of our analysis. As they are relevant to our method, we also present new calibrations of the dark matter halo mass function, halo mass accretion histories, and halo-subhalo merger rates out to z = 8. We also provide new compilations of cosmic and specific star formation rates; more recent measurements are now consistent with the buildup of the cosmic stellar mass density at all redshifts. Implications of our work include: halos near 10 12 M /circledot are the most efficient at forming stars at all redshifts, the baryon conversion efficiency of massive halos drops markedly after z ∼ 2 . 5 (consistent with theories of cold-mode accretion), the ICL for massive galaxies is expected to be significant out to at least z ∼ 1 -1 . 5, and dwarf galaxies at low redshifts have higher stellar mass to halo mass ratios than previous expectations and form later than in most theoretical models. Finally, we provide new fitting formulae for star formation histories that are more accurate than the standard declining tau model. Our approach places a wide variety of observations relating to the star formation history of galaxies into a self-consistent framework based on the modern understanding of structure formation in Λ CDM. Constraints on the stellar mass-halo mass relationship and star formation rates are available for download online.  \nSubject headings: dark matter - galaxies: abundances - galaxies: evolution - methods: N-body simulations', '1. INTRODUCTION': "Constraining the buildup of stellar mass in galaxies provides fundamental constraints on galaxy formation models. An ever growing number of galaxy observations have been made, covering a time period from 500 Myr after the Big Bang (e.g. Zheng et al. 2012) to the present day. A model of galaxy formation that can match observations over this entire stretch of time would represent a significant aid to our understanding. So far, matching the evolution of the stellar masses and star formation rates of galaxies over this entire epoch with has proved difficult (see, e.g., Lu et al. 2012; Borgani & Kravtsov 2011; Weinmann et al. 2012). Despite challenges of reproducing these observations with detailed models of galaxy formation, significant progress has been made in recent years with empirical models that connect the evolution of galaxy properties to the evolution of dark matter halos.  \nOver the past decade, a range of studies have associated galaxies with dark matter halos at a given epoch, using a variety of techniques, including Halo Occupation Distribution modeling (e.g., Berlind & Weinberg 2002; Bullock et al. 2002) the Conditional Luminosity Function modeling (e.g., Yang et al. 2003), and variants of the abundance matching technique (e.g., Colín et al. 1999; Kravtsov & Klypin 1999; Neyrinck et al. 2004; Kravtsov et al. 2004; Vale & Ostriker 2004, 2006; Tasitsiomi et al. 2004; Conroy et al. 2006; Shankar et al. 2006; Berrier et al. 2006; Marín et al. 2008; Guo et al. 2010; Moster et al. 2010; Behroozi et al. 2010). The simplest abundance matching models assign the most  \nmassive observed galaxies in rank order to the largest halos in an equal simulation volume. With only slight modifications (e.g., using the peak mass for satellite halos), this technique accurately reproduces the redshift- and scale-dependent clustering of galaxies (e.g., Conroy et al. 2006; Reddick et al. 2012).  \nBecause the cosmological model provides a prediction for the buildup of dark matter halos, one can combine knowledge of the galaxy-halo connection at different epochs with knowledge of the mass accretion and merger histories of halos to constrain the buildup of stellar mass in galaxies over cosmic time. This was first done in a comprehensive way by Conroy & Wechsler (2009), who provided an empirical constraint on the star formation histories and stellar mass growth in galaxies from z = 2 to the present. This approach has been explored in several studies using variants of the techniques (Zheng et al. 2007b; White et al. 2007; Firmani & Avila-Reese 2010; Leitner 2012; Béthermin et al. 2012; Wang et al. 2013; Moster et al. 2013), which together provide important constraints on the basic picture for the buildup of stars in galaxies and their connection to dark matter halos.  \nThese studies represented important advances, albeit with shortcomings. For example, most of these studies only modeled galaxy evolution from z ∼ 2 to the present due to perceived conflicts between the integrated cosmic star formation and the cosmic stellar mass density at z > 1 (Wilkins et al. 2008; Hopkins & Beacom 2006). In addition, because passive stellar mass loss (e.g., supernovae of massive stars) de-  \nnds on galaxies' star formation histories, these studies had to make assumptions about the ages of stars already formed at z > 2. Finally, these studies paid limited attention to the vast array of observational uncertainties as well as modeling uncertainties affecting their derived constraints; indeed, many results are presented without error bars (Conroy & Wechsler 2009; Firmani & Avila-Reese 2010; Béthermin et al. 2012; Wang et al. 2013).  \nHere we present a new technique and a new compilation of observations, aimed squarely at resolving these issues. On the observational side, new observations of galaxies at very high redshifts (e.g. Bouwens et al. 2011; Bradley et al. 2012; Bouwens et al. 2012; McLure et al. 2011) have made it possible to place constraints on galaxy formation all the way to z ∼ 8. At lower redshifts, we present the first constraints from stellar mass functions based on the PRIMUS survey, as well as from a new analysis at z = 0 based on SDSS and GALEX data (Moustakas et al. 2013). These are obtained from a consistent methodology from z = 0 to z = 1, and thus alleviate many concerns about matching inconsistently-identified galaxies at different epochs. Using these combined data sets, we find that observations of the evolution of galaxy star formation rates and stellar masses can now be reconciled (see also Bernardi et al. 2010 and Moster et al. 2013).  \nOur new method constrains the galaxy-halo relation using observed galaxy stellar mass functions as well as specific star formation rates and cosmic star formation rates. As with previous studies, we match observed galaxies to halos, but the additional information on star formation rates allows us to break degeneracies and directly constrain the buildup of the intracluster light, as well as the amount of stars that can transfer from satellites to the central galaxy during mergers. We account for a number of statistical and systematic effects, including uncertainties from stellar population synthesis models, dust models, star formation history models, the faint-end slope of the stellar mass function, observational completeness, scatter between stellar mass and halo mass.  \nGiven a parametrization for the intrinsic stellar mass-halo mass relationship as well as for these observational effects, we use a Markov Chain Monte Carlo (MCMC) method to find the allowed posterior distribution. This is important not only for addressing concerns that our results may be biased by limited observational constraints at high redshifts, but it also gives us the power to determine which observations would best improve the resulting constraints on the relationship between stellar mass, star formation, and halo mass. The direct results of this analysis are constraints on the distribution of galaxy stellar masses as a function of halo mass and redshift. From these, many other constraints relevant to galaxy formation are derived, including the average star formation rate as a function of halo mass and redshift, the average star formation history in galaxies at a given epoch as a function of galaxy stellar mass or host halo mass, the instantaneous baryon conversion efficiency of galaxies as a function of mass and redshift, the buildup of the intracluster light, and the evolution of the stellar mass to halo mass ratio for progenitors of today's galaxies, all including full uncertainties and covering a redshift range from z = 0 to z = 8.  \nWe provide a broad overview of our methodology in §2, including our parametrization of the stellar mass - halo mass relation and the relevant uncertainties, with additional details in Appendices A-B. We discuss the observational data sets relevant to our method in §3, with special attention to new measurements of cosmic star formation rates. We discuss  \nthe simulations that we use in §4, along with recalibrations of the halo mass function (Appendix G), halo mass accretion rates (Appendix H), and subhalo merger/disruption rates (Appendix I). We present our main results in §5, with discussion in §6 and a summary of our conclusions in §7.  \nThroughout this work, we assume a Chabrier (2003) initial mass function, the stellar population synthesis model of Bruzual & Charlot (2003), and the dust model employed in Blanton & Roweis (2007). We convert data sets from other papers to these models as necessary. We additionally assume a flat, Λ CDM cosmology with parameters Ω M = 0 . 27, Ω Λ = 0 . 73, h = 0 . 7, ns = 0 . 95, and σ 8 = 0 . 82.", '2.1. Overview': 'Much previous work has gone into determining the relation between halo mass and stellar mass as a function of redshift (e.g., Moster et al. 2013; Yang et al. 2012; Leauthaud et al. 2012; Behroozi et al. 2010; Moster et al. 2010; Guo et al. 2010; Conroy & Wechsler 2009; Wang & Jing 2010; Zheng et al. 2007b; Mandelbaum et al. 2006; Gavazzi et al. 2007; Yang et al. 2009a; Hansen et al. 2009; Lin & Mohr 2004). The redshift evolution of this relation is due to three contributing physical effects for stars in galaxies (new star formation, merging satellite galaxies, and stellar mass loss) and one effect for halos (continuing mass accretion). Merging galaxies and mass accretion are well-constrained by dark matter simulations, and stellar mass loss is well-constrained with the assumption of an initial mass function (IMF) for stars. The most uncertain of these effects, in relation to dark matter halos, are new star formation and stars ejected during galaxy mergers. With observational constraints on the redshift evolution of the stellar mass and the star formation rate as well as computed constraints on satellite mass loss, we may effectively constrain both of these factors.  \nSchematically, the star formation rate for a given galaxy over a given timestep is constrained by:', 'SFR · ∆ t =Expected Stellar Mass Now': "- -Remaining Stellar Mass from Previous Timestep\n- -(Stellar Mass from Mergers)  \n× (1 -Fraction Ejected) . (1)  \nIn this way, we can combine a specific assignment of stellar mass to halos with the evolution of the halo mass function over time, to derive star formation rates for halo trajectories. Our approach is to flexibly parametrize the possibilities for the stellar mass - halo mass relation, M ∗ ( Mh , z ), as well as the uncertainties affecting the remaining terms in Eq. 1. Then, using a Markov Chain Monte Carlo method, we can determine the allowed parameter space for M ∗ ( Mh , z ) by comparing the implied SFRs and stellar mass abundances to observations across a wide range of redshifts.  \nParametrizing M ∗ ( Mh , z ) separately for individual halos is beyond the scope of this work (although it will be explored in future papers). Here, we parametrize the median value of M ∗ for halos at a given mass and redshift (i.e., a M ∗ ( Mh , z )) as well as the scatter in stellar mass at fixed halo mass as a function of redshift. By taking a summation of Eq. 1 over halos in a specific mass and redshift bin, we can derive the total amount of stars formed in that bin. Thus, dividing by the number of halos in that bin gives the average star formation rate, which is similar to approaches taken in previous  \nworks (Conroy & Wechsler 2009; Leitner 2012; Wang et al. 2013; Moster et al. 2013). This approach is accurate except in terms of small second-order effects that come from different stellar mass histories having a range of stellar mass loss rates. This aspect is discussed in detail in Appendix A. Further discussion of how individual galaxy star formation histories can differ from the average is presented in §5.8.  \nIn this section, we present our method for parametrizing the stellar mass - halo mass relation (§2.2.1), our parametrization of the uncertainties affecting observational data (§2.2.2), our model for stellar mass loss (§2.3), and our model for ejection of stars into the ICL (§2.4). To streamline the presentation of this section, we present the methodology for incorporating the effects of stellar mass accreted in mergers (which involves calibrating merger rates from simulations to z = 8) and for calculating observables (such as the stellar mass function and specific star formation rates) to Appendices B and C. We present a summary of the methodology and details of the Markov Chain Monte Carlo method in §2.5.  \n2.2. Determining the Stellar Mass - Halo Mass Relation to z = 8  \n2.2.1. Intrinsic Relation  \nAs in Behroozi et al. (2010), we parametrize the stellar mass - halo mass (SMHM) relation at a given epoch. The most commonly-used function is a double power law (Yang et al. 2012; Moster et al. 2013), which has a characteristic halo mass ( M 1) and stellar mass ( M ∗ , 0), a low-mass slope ( α ), and a high-mass slope ( β ). Because we include data down to low stellar masses ( M ∗ ∼ 10 7 . 25 M /circledot ) and have tight observational constraints from the SDSS, we find that a double power-law cannot accurately fit the unique shape of the stellar mass function (SMF). As shown in Appendix D, a double power-law results in a stellar mass function off by as much as 0.1 dex at z = 0. For this reason, we choose a different, fiveparameter form as the fitting function for the SMHM relation. The chosen form retains the low-mass power-law slope ( α ) as well as a characteristic stellar ( M ∗ , 0) and halo mass ( M 1). The high-mass behavior is trickier to fit. We find the best match with a subpower (superlogarithm) function with index γ :  \nlog 10 ( M ∗ ( Mh →∞ )) ∝ ( log 10 ( Mh M 1 )) γ . (2)  \nIf γ is 1, this equation becomes an ordinary power law. As γ approaches 0, this equation approaches a logarithm. In between, this equation grows more slowly than any power law, but faster than any logarithm. While this equation can become multiply-valued for Mh < M 1, it is straightforward to mitigate this behavior as well as allowing a smooth connection to a power-law for Mh < M 1. The specific parametrization we adopt is: 1  \nlog 10 ( M ∗ ( Mh )) = log 10 ( /epsilon1 M 1) + f ( log 10 ( Mh M 1 )) -f (0) f ( x ) = -log 10 (10 α x + 1) + δ (log 10 (1 + exp( x ))) γ 1 + exp(10 -x ) . (3)  \nThis is a power law with slope -α for Mh /lessmuch M 1 and a subpower law with index γ for Mh /greatermuch M 1. The characteristic stellar mass to halo mass ratio is /epsilon1 at the characteristic halo mass  \nM 1. The maximum errors of this fit to the SMHM relation are about 0.025 dex, or roughly four times better than a double power-law fit, as discussed in Appendix D.  \nWe use the virial mass (as defined in Bryan & Norman 1998) to define the halo mass of central galaxies, and for satellites, we use the peak progenitor virial mass ( M peak). As this paper was being prepared, Reddick et al. (2012) found that the peak v max for halos is an even better proxy for stellar mass, at least at z = 0 . 05, in agreement with physical expectations that the depth of the halo potential well before it is impacted by stripping should be most correlated with the galaxy mass. The use of M peak is at present easier both conceptually (to consider mass ratios between galaxies and halos) and operationally (because it allows the use of several previously calibrated relationships in the literature), so we stick with this choice in the present work for simplicity. The use of M peak may result in slight underestimates of the clustering and the satellite fraction (see discussion in Reddick et al. 2012 for a comparison of these quantities in the local universe). However, §5.6 demonstrates that the differences between our results and the best-fitting stellar mass-halo mass relation of Reddick et al. (2012) are small, and we expect that they are within the current systematic errors.  \nThe redshift scaling of the relation poses a unique challenge. At z = 0, the parameter fits must be flexible enough to match the tight constraints from SDSS observations. At higher redshifts, the fit must be flexible enough to allow exploration of possible star formation histories allowed by higher uncertainties, but not so flexible that it allows significant over-fitting of the systematic and statistical errors in the data. Moreover, the scaling must be physical: obviously, negative star formation rates cannot be allowed, nor can stellar mass to halo mass ratios that exceed the halo's baryon fraction. These physical constraints help ensure basic sanity for the functional form where no observational data exists (for redshifts z > 8, as well as for faint galaxies).  \nBy necessity, the error bars we obtain will be sensitive to the choice of the redshift fit. This is not only so for the variation allowed in the stellar mass histories, but it is also true for the variation allowed in the nuisance parameters for systematic errors. Because the latter are often poorly constrained by available data, the flexibility of the redshift fit also influences how much of the available nuisance parameter space can be explored. For this reason, we have tried a large number of different redshift scalings, including cubic spline interpolation with many control points (which results in over-fitting) and many different choices of scaling parameters with a or z -many of which have unphysical behavior at very high redshifts. We have also tried fitting directly to abundance matching results (Appendix E), but no robust trends with redshift emerged.  \nThe choice that we settled on has several conventions. For many variables in Eq. 3, we have three parameters: one each for the low-redshift value, one for the scaling at intermediate redshifts (0 . 5 /lessorsimilar z /lessorsimilar 2), and one for scaling at high redshifts ( z /greaterorsimilar 2). We also apply an exponential shutoff to the redshift scaling towards lower redshift to isolate the highly constrained low-redshift parameters from higher-redshift ones. Our final parametrization is the following:  \nν ( a ) = exp( -4 a 2 ) log 10 ( M 1) = M 1 , 0 + ( M 1 , a ( a -1) + M 1 , zz ) ν  \nlog 10 ( /epsilon1 ) = /epsilon1 0 + ( /epsilon1 a ( a -1) + /epsilon1 zz ) ν + /epsilon1 a , 2( a -1) α = α 0 + ( α a ( a -1)) ν δ = δ 0 + ( δ a ( a -1) + δ zz ) ν γ = γ 0 + ( γ a ( a -1) + γ zz ) ν. (4)  \nAdditionally, there will be scatter in stellar mass at fixed halo mass ( ξ ). While current studies have indicated that this is approximately 0.16-0.2 dex at z = 0 (More et al. 2009; Yang et al. 2009b; Reddick et al. 2012) with no evidence for a mass trend at least down to 10 12 M /circledot in halo mass (Reddick et al. 2012), little is known about the redshift evolution of this scatter. For this work, we parametrize the possible evolution of scatter with redshift via a two-parameter scaling:  \nξ = ξ 0 + ξ a ( a -1) , (5)  \nand we take the prior on ξ 0 to be 0 . 20 dex ± 0 . 03 dex, as determined by Reddick et al. (2012) using a combination of constraints from the correlation function and conditional stellar mass function.", '2.2.2. Observational Systematics': "A full discussion of the systematic uncertainties affecting stellar mass functions may be found in Behroozi et al. (2010). These include uncertainties from the Initial Mass Function, the stellar population synthesis model, the dust model, the star formation history model, sample variance, Eddington bias, redshift errors, and magnification bias. Most of these effects result in a constant systematic offset ( µ ) in stellar masses:  \nlog 10 ( M ∗ , meas M ∗ , true ) = µ. (6)  \nHowever, some effects can have different behavior as a function of stellar mass, largely on account of assumptions for passive galaxies (e.g., declining star formation rates and low dust fractions) not being true for active galaxies. We account for these effects by introducing an additional parameter ( κ ) to capture offsets in the stellar masses of active galaxies:  \nlog 10 ( M ∗ , meas , active M ∗ , true , active ) = µ + κ. (7)  \nFor the fraction of galaxies that are quiescent, we use a fitting formula for recent measurements from Brammer et al. (2011), corrected to the stellar mass estimates we use in this paper:  \nf passive( M ∗ , meas , z ) = [ ( M ∗ , meas 10 10 . 2 + 0 . 5 z M /circledot ) -1 . 3 + 1 ] -1 . (8)  \nAt low redshifts, most galaxies below 10 10 . 2 M /circledot in stellar mass are active; beyond z = 3, nearly all galaxies are considered active.  \nAs in Behroozi et al. (2010), we do not model uncertainties in the IMF in this work. The most significant remaining uncertainties come from the SPS model ( ∼ 0 . 1dex), the dust model ( ∼ 0 . 1dex), and the star formation history ( ∼ 0 . 2dex). By definition, passive galaxies have had no recent star formation, so that the star formation history uncertainties are substantially smaller. For that reason, we take the priors on µ to be a log-normal distribution centered at zero with width 0.14 dex (combined SPS and dust model errors), but the priors on κ to be of width 0.24 dex (all three sources combined). The  \nfunctional form of how µ and κ evolve with redshift is unknown. We therefore adopt the following fiducial formulae:  \nµ = µ 0 + ( a -1) µ a (9)  \nκ = κ 0 + ( a -1) κ a . (10)  \nAt higher redshifts ( z ∼ 2 . 3), Muzzin et al. (2009) suggest that uncertainties from SPS modeling become more important, on the order of 0.2 dex instead of 0.1 dex. We thus set the widths of the priors on µ a and κ a to be 0.22 dex and 0.30 dex, respectively, which take these increased uncertainties into account.  \nBecause stellar populations cannot be fully constrained with limited photometric information, stellar mass estimates for individual galaxies have intrinsic scatter relative to the true galaxy stellar mass. This causes an Eddington bias (Eddington 1913, 1940) in the stellar mass function: on the high-mass end of the SMF, many low-mass galaxies are upscattered, but there are a limited number of higher-mass galaxies that can be downscattered. This results in a net increase in the observed numbers of high stellar-mass galaxies. This effect is best estimated for each redshift range and modeling technique individually; while few authors have corrected for this error, many at least provide an estimate of its effect. As in Behroozi et al. (2010), we find that these estimates depend most significantly on the redshift. We thus model the distribution in the observed stellar mass estimates (compared to the true stellar masses) as a log-normal Gaussian with mean 0 and a redshift-dependent standard deviation given by:  \nσ ( z ) = σ 0 + σ zz . (11)  \nFollowing Conroy et al. (2009); Behroozi et al. (2010), we take σ 0 = 0 . 07 and we take the prior on σ z to be σ z = 0 . 04 ± 0 . 015, consistent with estimates from the literature (Conroy et al. 2009; Kajisawa et al. 2009; Pérez-González et al. 2008; Marchesini et al. 2009; Stark et al. 2009; Caputi et al. 2011; Lee et al. 2011; Marchesini et al. 2010).  \nFinally, there are some concerns that, due to the Lyman break techniques used to detect high-redshift galaxies, not all high-redshift galaxies may be detected. If galaxies have bursty star formation, or if some fraction of star forming galaxies have extreme quantities of dust, then not all galaxies will be detected in Lyman break surveys (see also Stark et al. 2009; Lee et al. 2009). We parametrize the stellar mass completeness of high-redshift surveys with a two-parameter fit, governing the amplitude ( A ) and redshift onset ( zc ) of incompleteness:  \nci ( z ) = 1 -A exp( zc -z ) + 1 . (12)  \nBecause there is little reason to believe that surveys are missing a significant fraction of galaxies at z < 1, we require zc > 0 . 8, and we set the completeness fraction to 1 for z < 1. This amounts to adopting the following formula for the completeness:  \nc ( z ) = { 1 if z < 1 ci ( z ) + (1 -ci (1)) if z > 1 (13)  \nWe parametrize the effect of this incompleteness on the observed SFRs by a single parameter b , which sets the fraction of 'bursty' vs. 'dusty' star formation. Fully bursty star formation ( b = 1) would have no effect on the total observed cosmic SFR, but the observed SSFRs would be boosted by the  \nincompleteness factor compared to the SSFR averaged over longer periods of time. Fully dusty star formation ( b = 0) would lower the observed cosmic SFR by the same factor as the cosmic SM density. Presumably, a moderate amount of unobscured star formation would be enough to render a high-redshift galaxy observable-hence, dustiness will impact galaxy completeness only if star formation is largely obscured. Given that models for SFRs from observed galaxies already attempt to correct for dust obscuration, the main effect of dustiness on SSFRs is that the SSFRs for very dusty galaxies will not be observed. Therefore, in our current model, we assume that this has no effect on the average SSFR.", '2.3. Star Formation Histories': 'Under our assumption of a Chabrier (2003) IMF and Bruzual & Charlot (2003) stellar evolution tracks, we use the FSPS package (Conroy et al. 2010; Conroy & Gunn 2010) to calculate the rate of stellar mass loss. We find the fraction of mass lost from a single stellar population as a function of the time since its formation to be well-fit by the following formula: 2  \nf loss( t ) = 0 . 05ln ( 1 + t 1 . 4 Myr ) (14)', '2.4. Extragalactic Light and Merging Galaxies': "From Eq. 1, galaxies can build up stars either through mergers or internal star formation. Schematically, this implies that the unknown amount of stars deposited in mergers can be derived given separate constraints on the stellar mass growth of galaxies (i.e., from stellar mass functions) and their star formation rates:  \n∆ SM deposited = ∆ SM central -SFR central ∆ t , (15)  \nwhere ∆ SM deposited is the amount of stellar mass deposited from merging satellites, ∆ SM central is the stellar mass growth of the central galaxy, and SFR central is the star formation rate of the central galaxy over the chosen period ∆ t . (Note that the real equation is slightly more complicated due to stellar mass loss from passive evolution; see Appendix B).  \nWhen a galaxy merger occurs, the stars associated with the satellite galaxy may either be deposited onto the central galaxy or be ejected into the intrahalo light (IHL), also called the intracluster light (ICL) for galaxy clusters. Some care is necessary when using these terms, because there exist two separate definitions of the ICL/IHL from simulations (i.e., stars not bound to the main galaxy) and from observations (stars not counted in the light profile of the main galaxy). In this work, we use 'intracluster light' (and 'ICL') to mean any remnant stars from past galaxy mergers which are not counted as part of the light profile from the main galaxy, regardless of their boundedness. We do not make a distinction between IHL and ICL; we use one term ('ICL') regardless of the size of the host halo, as this usage is by far more common in the literature.  \nBecause we can constrain the total amount of stars in merging satellites through the halo merger rate (Appendix I) and the stellar mass - halo mass relation, we can also derive the amount of stellar mass deposited into the ICL ( ∆ ICL ):  \n∆ ICL = ∆ SM incoming -∆ SM deposited , (16)  \nwhere ∆ SM incoming is the total amount of stars in satellites that merge within ∆ t . Note that this quantity includes all stars ever ejected from merging galaxies. Some fraction of these ejected stars may have speeds higher than the escape velocity of the surrounding dark matter halo and will be scattered to very large distances; however, this fraction is expected to be on the order of a few percent or less (Behroozi et al. 2012).  \nThis simple picture is made more complicated by inconsistencies in how galaxies and the ICL are separated. The total luminosity recovered and the fraction of it attributed to the galaxy are both dependent on surface-brightness limits, sky subtraction methods, and galaxy fitting methods, which are different for the different data sets we use and for different redshift ranges in the data sets themselves. Over most mass ranges and redshifts, this is not a problem: the rate of star formation in many galaxies results in an overall change in stellar mass between redshifts that is much larger than that attributable to a change in the luminosity modeling method. However, for massive galaxies ( M ∗ > 10 11 M /circledot ) at low redshifts ( z < 1), this is a more significant issue. For example, including 5% more light between z = 0 . 3 and z = 0 . 1 in a massive galaxy (0.05 magnitudes) would result in an equal buildup in stellar mass as would a specific star formation rate of about 2.5 × 10 -11 yr -1 (the expected SSFR for such galaxies; Salim et al. 2007). As the true SSFR drops to 10 -12 yr -1 , the galaxy models must be self-consistent to about 0.2% (0.002 magnitudes) to avoid introducing a comparable error, which is well beyond current calibration methods (Bernardi et al. 2010).  \nThis said, changes in the mass of the galaxy because of definitional issues are completely degenerate with incoming mass from mergers in Eqs. 15 and 16: ∆ SM deposited is replaced by the sum of the definitional and deposited stellar mass changes ( ∆ SM definitional + ∆ SM deposited) everywhere it appears. Notably, Eq. 16 still gives us a way to robustly determine the buildup of stars in the ICL.  \nIf we used Eqs. 16 and 15 directly, the star formation rate would have to be parametrized separately from the stellar mass growth. Instead, we rewrite Eq. 15 in terms of the fraction of stellar mass growth that comes from the star formation rate:  \n∆ SM deposited + ∆ SM definitional = (1 -fSFR ) ∆ SM central (17)  \nfSFR = SFR central ∆ t ∆ SM central . (18)  \nThis is useful because there are strong physical priors on what fSFR can be. Low-mass galaxies do not have enough incoming stellar mass in satellites to account for a significant part of their stellar mass growth, so fSFR asymptotes to 100% for such galaxies. On the other hand, the highest-mass galaxies experience almost no internal star formation (Salim et al. 2007), meaning that fSFR must approach 0 for such galaxies. So, we can approximate the halo mass dependence of fSFR using a double power law as  \nfSFR ( Mh ) = [ ( Mh Mh , ICL ) β + 1 ] -1 . (19)  \nWe allow Mh , ICL to be redshift-dependent using a twoparameter fit:  \nMh , ICL = Mh , ICL , 0 + ( a -1) Mh , ICL , a . (20)  \nRather than add extra parameters for β , we make the assumption that fSFR (10 16 M /circledot ) is 1%. This is equivalent to the constraint  \nβ = log 10 (99) 16 -log 10 ( Mh , ICL ) . (21)", '2.5. Methodology Summary': 'Although the equations involved are somewhat complicated, the logical steps involved in our method are straightforward and are shown visually in Fig. 1. These steps include:  \n- 1. Choose parameters for the stellar mass - halo mass relation as well as observational and methodology uncertainties (see Table 1 for full list; see Table 2 for adopted priors) using a Markov Chain Monte Carlo method.\n- 2. Use the chosen stellar mass - halo mass relation to populate halos with galaxies, and use merger rates and mass accretion histories calculated from dark matter simulations to infer galaxy growth rates.\n- 3. From the previous step, calculate the average star formation rates as a function of halo mass and redshift, and use those to calculate the average specific star formation rate and the cosmic star formation rate. In addition, using the abundances of halos, calculate the stellar mass function.\n- 4. Apply corrections for observational errors and biases to the derived stellar mass function and star formation rates.\n- 5. Compare with observational measures of the stellar mass function and star formation rates (i.e., sum χ 2 errors for each data point) to calculate the likelihood that the chosen stellar mass - halo mass relation matches observed results (i.e., exp( -0 . 5 χ 2 )).\n- 6. Return to step #1 until the MCMC algorithm has converged.  \nTo ensure convergence on a space with such a large number of parameters, we use the Adaptive Metropolis exploration method (Haario et al. 2001). Although the algorithm of Dunkley et al. (2005) indicates that we converge after only 5 × 10 5 points, we continue to run for 4 × 10 6 total points to minimize the chance of unexplored regions and to ensure that the adaptive updates of the step covariance matrix converge and do not bias our final results.', '3. OBSERVATIONAL CONSTRAINTS AND UNCERTAINTIES': "We use three types of observational constraints on star formation in galaxies: the stellar mass function (SMF), the specific star formation rate of galaxies (SSFR, star formation rate per unit stellar mass), and the cosmic star formation rate (CSFR). Stellar mass functions are discussed in §3.1, cosmic star formation rates are discussed in §3.3, and specific star formation rates are discussed in §3.2.  \nTo ensure a fair comparison between the different data sets, we convert each measurement to the same set of systematic assumptions. For the initial mass function (IMF), we convert all results to that of Chabrier (2003). We do not consider evolving IMFs, nor do we consider the uncertainties for a nonuniversal IMF, as both uncertainties are beyond the scope of  \nthis paper. For the stellar population synthesis and dust models, we convert all results to the models of Bruzual & Charlot (2003) and Blanton & Roweis (2007). In all cases, we convert authors' supplied units to physical units under the assumption that h = 0 . 7.", '3.1. Stellar Mass Functions': 'We combine many overlapping data sets (Baldry et al. 2008; Moustakas et al. 2013; Pérez-González et al. 2008; Mortlock et al. 2011; Marchesini et al. 2009; Stark et al. 2009; Bouwens et al. 2011) to constrain the evolution of the stellar mass function from z = 0 to z = 8. These are shown in Table 3.  \nDifferent observational challenges present themselves in each redshift range. At low redshifts, the statistics are sufficient that systematic errors in calculating stellar masses are the largest source of uncertainty (Behroozi et al. 2010). At higher redshifts, this is no longer the case; statistical and sample variance errors can be of equal magnitude to estimated systematic errors. Nonetheless, under the assumption of a smoothly varying stellar mass - halo mass (SMHM) relation, our approach (which simultaneously constrains the SMHM relation against many stellar mass functions at different redshifts) can constrain the statistical errors better than the individual error bars on stellar mass functions might suggest. We refer the interested reader to Behroozi et al. (2010) for a discussion of the most common uncertainties affecting observation constraints on the stellar mass function and to §2.2.1 for our parametrization of their effects. In this section, therefore, we limit ourselves to discussing special concerns relevant to using stellar mass functions over this broad range of redshifts.  \nQuestions about the accuracy of stellar mass estimates over this redshift range come from several previously published results that suggest that the evolution in the cosmic SFR density is inconsistent with the estimated cosmic stellar mass density at z > 1 (Nagamine et al. 2006; Hopkins & Beacom 2006; Pérez-González et al. 2008; Wilkins et al. 2008). One possible explanation is an evolving or a non-universal IMF; a number of different lines of evidence exist in support of this theory (e.g., van Dokkum & Conroy 2010; Lucatello et al. 2005; Tumlinson 2007a,b; van Dokkum 2008). Another possible reason is that star formation at high redshift could be bursty or dust-obscured, yielding an incomplete census of galaxies in Lyman-break surveys (Lee et al. 2009; Stark et al. 2009). On the other hand, Reddy & Steidel (2009) offer a simpler explanation. They appeal to luminosity-dependent reddening corrections in the ultraviolet luminosity functions at high redshift as well as steeper faint-end slopes for stellar mass functions, and demonstrate that the purported discrepancy then largely vanishes.  \nIn our analysis of recent cosmic SFR literature (§3.2), we find that newer estimates of the cosmic SFR are systematically lower at z > 3 than the data in Hopkins & Beacom (2006) (see Fig. 2). A similar conclusion is reached in Bernardi et al. (2010) and Moster et al. (2013). In addition, deep probes of the stellar mass function at lower redshifts ( z < 1 . 5) have resulted in findings of steep faint-end slopes (Baldry et al. 2008; Drory et al. 2009; Mortlock et al. 2011). We find this to be strong evidence in support of the explanation in Reddy & Steidel (2009). Nonetheless, we also consider simple models for the effects of bursty or dusty modes of star formation; an analysis of models for an evolving or non-universal IMF is beyond the scope of this paper.  \nThe highest redshifts present unique challenges. Herein,  \nFIG. 1.- Visual summary of the methodology used to constrain the stellar mass - halo mass relation.  \n<!-- image -->  \nTABLE 1 TABLE OF PARAMETERS  \nfor redshifts 7 < z < 8 . 5, we convert the UV luminosity function given in Bouwens et al. (2011) to a stellar mass function according to the recipe in González et al. (2011). It has been suggested that this conversion is reasonable given the blue UV continuum slopes of such high-redshift galaxies, and follow-up observations from IRAC that imply that these early galaxies are nearly dust-free (Labbé et al. 2010). At the same time, concerns arise because the UV luminosity function presented in Bouwens et al. (2011) is steep ( α ≈ -1 . 6 to -2 . 0, instead of the shallower α ≈ -1 . 5 to -1 . 7 at lower redshifts). One interpretation of this result is that assumptions for galaxy morphology (and, therefore, the ability to measure estimate the full UV flux or half-light radii) can bias the calculated UV luminosity function to have a steeper faintend slope (Grazian et al. 2011). Indeed, any systematic effect that interferes more with brighter galaxies than with dimmer ones (e.g., dust) will tend to cause an overestimation of the slope of the UV luminosity function. González et al. (2011) assume that luminosity-dependent corrections for the UV-tostellar-mass conversion do not evolve substantially from z = 4 to z = 7. This may in fact overestimate the dust corrections,  \nan error that would lower the faint-end slope, partially correcting for a failure to account for galaxy morphology. In addition, recent work has shown that many estimates of stellar masses at z > 5 may have been overestimated due to nebular line emission (Stark et al. 2013); however, the conversions in González et al. (2011) do not share this problem (D. Stark, priv. comm.). As the use of UV-only stellar masses therefore remains controversial, we present results both with and without the converted Bouwens et al. (2011) stellar mass functions for z = 7 -8 in Appendix J.', '3.2. Cosmic Star Formation Rates': "To assemble constraints on cosmic star formation rates, we conducted a comprehensive literature search on the astro-ph arXiv for papers posted within the past 6 years (2006-2012); details of the search are in Appendix F. The selected papers are summarized in Table 4. Many different types of surveys were conducted, including estimates of the cosmic SFR from narrowband (e.g., H α ), broadband (UV-IR), and radio (1.4 GHz) surveys.  \nAs shown in Fig. 2, newer estimates of the cosmic SFR are  \nTABLE 2 TABLE OF PRIORS  \nNOTE. G ( x , y ) denotes a Gaussian distribution with center x and width y . U ( x 1 , x 2) denotes a uniform distribution from x 1 to x 2. The remaining parameters have no explicit priors; M 1 and /epsilon1 are explored in logarithmic space, whereas α , δ , and γ are explored in linear space.  \nTABLE 3 OBSERVATIONAL CONSTRAINTS ON THE STELLAR MASS FUNCTION  \nNOTE. - Letters correspond to the following corrections made to the published results: I (Initial Mass Function), S (Stellar Population Synthesis model), D (Dust model), V (Sample Variance), B (Surface Brightness incompleteness), C (Cosmology, specifically h ), U (UV to stellar mass conversion according to González et al. 2011). The local results ( z < 0 . 2) in the Moustakas et al. (2013) mass functions are taken from the Sloan Digital Sky Survey and cover an area of 2505 deg 2 .  \nsystematically and substantially lower than previous compilations (such as Hopkins & Beacom 2006) at z > 3. At high redshifts, the predominant method of probing star formation is through UV emission; however, conversion to a star formation rate requires quantifying dust obscuration systematics. Previously, Hopkins & Beacom (2006) assumed fixed dust obscuration for z > 2 . 5. However, new measurements have shown that the amount of dust is likely to decrease substantially with redshift beyond z ∼ 3 (Bouwens et al. 2012; Reddy & Steidel 2009). The inferred star formation rates thus were likely to have been overcorrected in Hopkins & Beacom (2006), although substantial uncertainty about the cosmic star formation rate at high redshifts still remains.  \nBecause most survey authors do not fully consider systematic errors, most reported error estimates (1-10%) are too small to explain the observed variance in published estimates of the cosmic star formation rate (see Fig. 2). Therefore, instead of using authors' estimates for systematic uncertainties, we estimate the true systematic errors by computing the variance of our collection of published SFR estimates, as detailed in Appendix F. Our results are given in Table 7. We find the average systematic error to range from 0.13 dex at z = 0 to 0.27 dex for z > 3.", '3.3. Specific Star Formation Rates': "We conducted a literature survey in an almost identical manner as for the cosmic star formation rate to assemble a  \nFIG. 2.- Cosmic star formation rates reported in the past six years have been dramatically lower at high redshifts ( z > 2) than previous results. Top panel: our new compilation of cosmic star formation rates, which excludes those reported before 2006 (see Table 4), along with a best-fitting double power law. Bottom panel: the compilation of cosmic star formation rates reported in (Hopkins & Beacom 2006) along with a best-fitting double power law (black dashed line) compared to the best-fitting double power law from our new compilation (blue dashed line). The largest differences at high redshifts result from changing assumptions about the amount of dust present at z > 3. Fits (Eq. F1) to the new and older data sets are given in Table 6. Error bars are not shown since publications have often dramatically underestimated the true magnitude of the systematic errors (see §3.2).  \n<!-- image -->  \nz  \n<!-- image -->  \nz  \ncollection of specific star formation rate (SSFR) constraints as a function of stellar mass, again in October 2012 (see Appendix F). The selected papers are summarized in Table 5. Most authors do not consider the full range of systematics that could affect their stellar masses or SFR estimates. To account for this we take a similar approach as in §3.2 and estimate  \nTABLE 4 OBSERVATIONAL CONSTRAINTS ON THE COSMIC STAR FORMATION RATE  \nNOTE. - Letters correspond to the following corrections made to the published results: I (Initial Mass Function), A (Intrinsic Scatter Correction; see Appendix F). The technique of Le Borgne et al. (2009) (parametric derivation of the cosmic SFH from counts of IR-sub mm sources) uses multiple surveys with different areas.  \nTABLE 5 OBSERVATIONAL CONSTRAINTS ON THE SPECIFIC STAR FORMATION RATE  \nNOTE. - The letter 'I' corresponds to an Initial Mass Function correction; 'N' corresponds to a correction for nebular emission lines (Stark et al. 2013). Daddi et al. (2007) report using a fraction of the GOODS-N and GOODS-S fields, but did not report an exact area. Schaerer & de Barros (2010) use a combination of fields, the largest of which is COSMOS (with 2 deg 2 area). However, this has a comparatively bright limiting magnitude.  \nsystematic errors by calculating inter-publication variances. Details of our approach are given in Appendix F. We did not find evidence for mass or redshift trends in the variance, and have adopted a uniform systematic uncertainty of 0.28 dex for SSFRs.  \nBesides these observational constraints, we apply one modeling constraint to SFRs in high-mass halos at low redshifts. The slope of the stellar mass - halo mass relationship flattens  \ntowards high halo masses (Behroozi et al. 2010; Moster et al. 2010), meaning that central galaxy stellar mass becomes a poor indicator of the host halo mass. Consequently, observed SSFRs (reported as functions of stellar mass) offer poor constraints on the star formation rate of high-mass halos. Nonetheless, observations of galaxy clusters have indicated that the central galaxies generally have very low star formation rates (e.g. Donahue et al. 2010; Rawle et al. 2012). To account for this we apply a weak prior based on the Wetzel et al. (2012) group catalog. The recent star formation history of Mh > 10 14 M /circledot halos is averaged from z = 0 to z = 0 . 2 (corresponding to the range of the SDSS survey used in Wetzel et al. 2012); if it is above 1 M /circledot yr -1 , the the total χ 2 error of the fit is increased by ( 〈 SFH 〉 -1 M /circledot yr -1 ) 2 ( M /circledot yr -1 ) -2 .", '4. SIMULATION DATA': "We rely mainly on the high-resolution Bolshoi simulation, described in Klypin et al. (2011). Bolshoi follows 2048 3 ( ≈ 8 . 6 billion) particles in a comoving, periodic box with side length 250 h -1 Mpc from z = 80 to the present day. Its mass resolution (1 . 9 × 10 8 M /circledot ) and force resolution (1 h -1 kpc) make it ideal for studying the evolution of halos from 10 10 M /circledot (e.g., satellites of the Milky Way) to the largest clusters in the universe (10 15 M /circledot ). Bolshoi was run as a collisionless dark matter simulation with the Adaptive Refinement Tree Code (ART; Kravtsov et al. 1997; Kravtsov & Klypin 1999) assuming a flat, Λ CDM cosmology ( Ω M = 0 . 27, Ω Λ = 0 . 73, h = 0 . 7, σ 8 = 0 . 82, and ns = 0 . 95). These cosmological parameters are consistent with results from both WMAP5 (Komatsu et al. 2009) and the latest WMAP7+BAO+H0 results (Komatsu et al. 2011).  \nHalos in Bolshoi were identified using ROCKSTAR, a sevendimensional halo finder that uses phase space plus temporal information (Behroozi et al. 2013b). Merger trees were generated using a new algorithm presented in Behroozi et al. (2013c), which enforces physically consistent evolution of halo properties across timesteps to provide increased robustness against anomalies in the halo finder. Halo masses are calculated using spherical overdensities, according to the virial overdensity criterion of Bryan & Norman (1998).  \nTo supplement results from Bolshoi where there exist concerns about statistics (i.e., for high-mass halos) or for initial conditions (as Bolshoi uses Zel'dovich/1LPT initial conditions), we also make use of two larger simulations. MultiDark follows the same number of particles in a much larger volume (1000 h -1 Mpc) with 7 times worse force resolution and 64 times worse mass resolution, but using the identical simulation code and cosmology (Riebe et al. 2011). Consuelo follows 1400 3 particles in a somewhat larger volume (420 h -1 Mpc) than Bolshoi; it has eight times worse force resolution and fourteen times worse mass resolution (McBride et al, in preparation; see also Behroozi et al. 2013c). 3 It has a slightly different cosmology ( Ω m = 0 . 25, Ω Λ = 0 . 75, h = 0 . 7, ns = 1 . 0, and σ =0 . 8) and was run using the GADGET-2 code (Springel 2005). Notably, however, it uses the more accurate 2LPT initial conditions instead of Zel'dovich initial conditions, so it serves as a way to correct for inaccuracies in the high-redshift mass function for Bolshoi.  \nResults from these three simulations are used to calibrate the halo mass function (Appendix G), halo mass accretion histories (Appendix H), and subhalo disruption rates (Appendix  \nFIG. 4.- The best fitting model (red line) and posterior one-sigma distribution (red shaded region) for the evolution of the specific star formation rate from z = 0 to z = 8, compared to observational estimates (black points).  \n<!-- image -->  \nO  \n<!-- image -->  \nz  \nFIG. 3.Top panel: Evolution of the stellar mass function from z = 0 to z = 8 in the best fitting model (colored lines), compared to observations (points with error bars; for clarity not all data is shown). Bottom panel: Observational constraints on the cosmic star formation rate (black points), compared to the best-fit model (red solid line) and the posterior one-sigma distribution (red shaded region).  \nI) used in this work.", '5. RESULTS': 'The method presented above results in a posterior distribution for the set of parameters describing models that match observed stellar mass functions, specific star formation rates, and cosmic star formation rates from z = 0 to z = 8. All data results in this paper are available for download online. 4 Our best-fitting parameters with one-sigma limits are as follows:  \nIntrinsic Parameters:  \nlog  \n<!-- image -->  \nz  \nρ 0 . 5 = 0 . 799 + 0 . 028 -0 . 355  \nSystematic Parameters:  \nν =exp( -4 a 2 ) log 10 ( /epsilon1 ) = -1 . 777 + 0 . 133 -0 . 146 + ( -0 . 006 + 0 . 113 -0 . 361 ( a -1) + ( -0 . 000 + 0 . 003 -0 . 104 ) z ) ν + -0 . 119 + 0 . 061 -0 . 012 ( a -1) log 10 ( M 1) = 11 . 514 + 0 . 053 -0 . 009 + ( -1 . 793 + 0 . 315 -0 . 330 ( a -1) + ( -0 . 251 + 0 . 012 -0 . 125 ) z ) ν α = -1 . 412 + 0 . 020 -0 . 105 + (0 . 731 + 0 . 344 -0 . 296 ( a -1)) ν δ =3 . 508 + 0 . 087 -0 . 369 + (2 . 608 + 2 . 446 -1 . 261 ( a -1) +-0 . 043 + 0 . 958 0 . 071 z ) ν γ =0 . 316 + 0 . 076 -0 . 012 + (1 . 319 + 0 . 584 -0 . 505 ( a -1) + 0 . 279 + 0 . 256 -0 . 081 z ) ν 10 ( Mh , ICL ) = 12 . 515 + 0 . 050 -0 . 429 + ( -2 . 503 -0 . 202 -2 . 078 )( a -1)  \nµ = -0 . 020 + 0 . 168 -0 . 096 + 0 . 081 + 0 . 078 -0 . 036 ( a -1) κ =0 . 045 + 0 . 110 -0 . 051 + ( -0 . 155 + 0 . 133 -0 . 133 )( a -1) ξ =0 . 218 + 0 . 011 -0 . 033 +-0 . 023 + 0 . 052 -0 . 068 ( a -1) σ =0 . 070 + 0 . 061 + 0 . 017 -0 . 008 ( z -0 . 1) ci ( z ) = 0 . 273 + 0 . 103 -0 . 222 (1 + exp(1 . 077 + 3 . 502 -0 . 099 -z )) -1 b =0 . 823 + 0 . 043 -0 . 629  \nOur total χ 2 error for the best-fit model from all sources (observational and theoretical) is 245. For the number of observational data points we use (628), the nominal reduced χ 2 is 0.4. While the true number of degrees of freedom is not  \n<!-- image -->  \nFIG. 5.Left panel: Average star formation rates as a function of halo mass and redshift. The overlaid white lines show average mass accretion histories for halos as a function of redshift for comparison. The grey area shows halos that would have a mass of > 10 15 . 5 M /circledot at z = 0 and therefore are not expected to exist. Right panel: Star formation histories (SFH) as a function of present-day halo mass and redshift, for galaxies at z = 0. This figure shows the historical star formation rate for stars in the galaxy at the present day. Since the contribution of stars from merging galaxies is so low, this is equivalent to the star formation rate traced along the white mass accretion trajectories in the left panel.  \n<!-- image -->  \n<!-- image -->  \nFIG. 6.Left panel: Average star formation rates for the galaxies in halos at a given halo mass and redshift (lines). Shaded regions indicate the one-sigma posterior distribution. Right panel: Average star formation histories as a function of halo mass and redshift (lines). Shaded regions indicate the one-sigma posterior distribution. Histories for 10 15 M /circledot halos are not shown as they are very similar to those for 10 14 M /circledot halos.  \n<!-- image -->  \n<!-- image -->  \nFIG. 7.Left panel: Evolution of the derived stellar mass as a function of halo mass. In each case, the lines show the mean values for central galaxies. These relations also characterize the satellite galaxy population if the horizontal axis is interpreted as the halo mass at the time of accretion. Error bars include both systematic and statistical uncertainties, calculated for a fixed cosmological model (see §4 for details). Right panel: Evolution of the derived stellar mass fractions ( M ∗ / Mh ) as a function of halo mass.  \n<!-- image -->  \nFIG. 8.- Evolution of the derived stellar mass fractions ( M ∗ ( z ) / Mh ( z )) as a function of halo mass at the present day. More massive halos used to have a significantly larger fraction of mass in stars, but the peak star formation efficiency has remained relatively constant to the present day.  \n<!-- image -->  \neasily computed (note for example that covariance matrices are not available for most stellar mass functions in the literature), this suggests that our best fit is a reasonable match to the data. This is shown visually in Fig. 3, which shows the evolution of the stellar mass function and the posterior distribution for the observed cosmic star formation rate and in Fig. 4, which shows the posterior distribution of specific star formation rates.  \nAs discussed in §2, each model in the distribution contains complete information on the average stellar content and formation history of halos as a function of mass and redshift. Herein, we focus on a few interesting implications, leaving more extensive coverage to a follow-up paper. In §5.1, we present derived star formation rates and star formation histories as a function of halo mass and redshift. Next, in §5.2, we discuss constraints on the stellar mass to halo mass relation, both as a function of historical halo mass and as a function of halo mass at z = 0. In §5.3, we compare the trajectories of the stellar mass and intracluster light (ICL) buildup in galaxies. Then, in §5.4, we discuss instantaneous baryon conversion efficiencies and how they relate to integrated baryon conversion efficiencies as well as stellar ages and formation times in §5.5. In §5.6, we compare our main results to those obtained by previous studies. We discuss the main effects of the uncertainties we have modeled in §5.7; finally, in §5.8, we present new fitting formulae for individual galaxy star formation histories relevant to observers.', '5.1. Star Formation Rates and Histories': "We show derived star formation rates as a function of halo mass in the left panel of Fig. 5, and the corresponding star formation histories for galaxies at z = 0 in the right panel. As discussed in §5.3, the contribution from merging galaxies to the central galaxy is small, so star formation histories for galaxies at z = 0 trace the star formation rate as a function of their progenitor's halo mass. Also, for this reason, star formation histories for galaxies at a given redshift zg are nearly the same as for the the z = 0 stellar populations, except truncated at z = zg . 5  \nWe show similar plots with one-sigma uncertainties in Fig. 6. The left-hand panel demonstrates that the star formation rate at fixed halo mass has been monotonically decreasing since very early redshifts. This rate of decrease is different for different halo masses. At moderate to high redshifts ( z > 2), larger halo masses generically have larger average star formation rates. However, at lower redshifts, the highest mass halos ( Mh /greaterorsimilar 10 14 M /circledot ) become so inefficient that they have lower star formation rates than group-scale (10 13 M /circledot ) halos or Milky-Way sized (10 12 M /circledot ) halos.  \nFrom the perspective of individual galaxies, it is more illuminating to look at the star formation history in the right panel of Fig. 6. Because halos continually gain mass over time, and do so more rapidly at early redshifts, the star formation history for galaxies is not monotonically decreasing. Instead, it increases rapidly with time, approximately as a power law in time. Depending on present-day halo mass, the galaxy's star formation rate reaches a peak at a redshift between z = 0 . 5 to z = 2 . 5 (higher redshifts for higher halo masses) and then decreases until the present day. The rate of decrease depends again on the halo mass, with high halo masses shutting off more rapidly than lower halo masses. Cluster-scale ( Mh /greaterorsimilar 10 14 M /circledot ) halos form most of their stars rapidly, at early times, whereas galaxies in Magellanic Cloudscale halos (10 11 M /circledot ) form stars over an extended period of time (see also §5.4).", '5.2. The Stellar Mass - Halo Mass Relation': "We show constraints on stellar mass - halo mass (SMHM) relation from z = 0 to z = 8 in the left panel of Fig. 7 and on the stellar mass - halo mass ratio in the right panel. As seen in our previous work (Behroozi et al. 2010), there is a strong peak in the stellar mass to halo mass ratio at around 10 12 M /circledot to at least z ∼ 4 and a weaker peak still visible to z ∼ 8. While the location of the peak appears to move to higher masses with increasing redshift (consistent with Leauthaud et al. 2012), the abundance of massive halos is also falling off with increasing redshift.  \nIn terms of dwarf galaxies ( Mh ∼ 10 10 M /circledot ), we only have constraints from observations at z = 0. These have been the subject of recent interest due to the finding of higher-thanexpected stellar mass to halo mass ratios in dwarf galaxies around the Milky Way (Boylan-Kolchin et al. 2012). However, these expectations have been set largely by the assumption that the stellar mass to halo mass ratio remains a scalefree power law below 10 11 M /circledot . As seen in Fig. 7, the low mass power-law behavior is broken below 10 11 M /circledot , corresponding with an upturn in the stellar mass function below 10 8 . 5 M /circledot (Baldry et al. 2008) (this result has also been seen by Kravtsov, in prep). This underscores the danger of assuming that faint dwarfs obey the same physical scaling relations as Magellanic Cloud-scale galaxies; moreover, it also is a strong argument against fitting the stellar mass - halo mass relation with a double power law (see discussion in Appendix D).  \nConcerning the range of allowed SMHM relations, the observational systematics are large enough that our results are marginally consistent with an unchanging SMHM relation from z = 6 to z = 0. Nonetheless, the feature with strongest significance is a gradual decrease in stellar mass in the median 10 11 M /circledot halo from z = 0 to z = 2, followed by an increase again for redshifts z > 6; also potentially indicated is an increase in stellar mass in the median Mh > 10 13 M /circledot halo from z = 0 to z = 2. The best fit SMHM relations at z = 7 and z = 8  \nz  \n<!-- image -->  \n<!-- image -->  \nz  \nFIG. 9.Left panel: Stellar mass growth histories of galaxies of different masses. Lines shown the amount of stellar mass remaining in central galaxies at the present day that was in place (in any progenitor galaxy) at a given redshift for our best fit model. Right panel: Amount of stellar mass in the intracluster light (ICL) remaining at the present day in place at a given redshift. Note that plots for 10 15 M /circledot halos are not shown as they are nearly identical to those for 10 14 M /circledot halos. Shaded regions in both panels shown the one-sigma posterior distribution.FIG. 10.- The fraction of stellar mass growth in galaxies due to in situ star formation (as opposed to growth by galaxy-galaxy mergers) as a function of halo mass and redshift.  \n<!-- image -->  \nare significantly different than at lower reshifts, with more stellar mass per unit halo mass. However, concerns about the reliability of the stellar mass functions at those redshifts (see §3.1) urge caution in interpreting the physical meaning of this result.  \nA useful perspective on these results can be obtained by considering the historical stellar mass to halo mass ratio of halos, as shown in Fig. 8. Despite the large systematic uncertainties, it is clear that halos go through markedly different phases of star formation. This evolution is most apparent for massive halos, as observations have been able to probe the properties of the progenitor galaxies all the way to z =8. Specifically, high-redshift progenitors of today's brightest cluster galaxies ( Mh ∼ 10 14 M /circledot were relatively efficient in converting baryons to stars-comparable to the most efficient galaxies today. However, between redshifts 2 -3, their efficiencies peaked, and thereafter they began to form stars less rapidly than their host halos were accreting dark matter. At the present day, such galaxies have an integrated star formation efficiency that is two orders of magnitude less than at their peak. The picture is less clear for progenitors of lowermass galaxies because current observations cannot probe their progenitors as far back. Nonetheless, their apparent behavior  \nin Fig. 7 of rising to a peak efficiency and later falling is consistent with all available data.", '5.3. Stellar Mass and Intracluster Light Growth Histories': "Our model constrains the buildup of stars in the intracluster light (see definition in §2.4) purely from observational galaxy data and measurements of the halo-halo merger rate in simulations (see also Watson et al. 2012 for an alternate method). In our best-fitting model, only 5% of stellar mass in mergers for 10 14 M /circledot halos is allowed to be deposited onto the central galaxy since z = 1, and only 10% for 10 13 M /circledot halos. For Milky Way-sized and smaller halos, this number rises rapidly to 70-80%. Yet, due to the sharply decreasing stellar mass to halo mass ratio for lower-mass halos, most of the incoming stellar mass will be in (rare) major mergers. Central galaxies are therefore relatively uncontaminated by stars from smaller recently-merged satellites. At higher redshifts, however, larger fractions of the stellar mass in merging galaxies are allowed to be deposited onto the central galaxy.  \nIn Fig. 9, we show the amount of galaxies' present-day stellar mass and intracluster light (ICL)/halo stars that was in place at a given redshift. The left-hand panel (stellar mass) shows that almost all stars in the central galaxy in present-day cluster-scale halos were in place at z = 2. However, since that time, a large number of their satellites have been disrupted into the ICL. Thus, for cluster scale halos, the stellar mass in the ICL exceeds the stellar mass in the central galaxy by a factor of 4-5, consistent with observations (e.g., Gonzalez et al. 2005). We note that our model predicts what may seem to be a large ICL fraction for Milky Way-sized galaxies (10 12 M /circledot ). However, the Milky Way is a special case. It has not had a major merger for ∼ 10-11 Gyr (Hammer et al. 2007); however, as noted above, only major mergers can contribute substantially to the ICL. A major merger 11 Gyr ago would have, however, contributed less than 3% of the present-day stellar mass of the Milky Way into the ICL; allowing for passive stellar evolution, this would result in less than 2% of the luminosity of the Milky Way coming from the intrahalo light, in excellent agreement with observations (Purcell et al. 2007).  \nFinally, in Fig. 10, we show the inferred fraction of stellar mass growth in galaxies coming from in situ star formation (as opposed to galaxy-galaxy mergers) as a function of halo mass and redshift. At all redshifts greater than 1, the vast majority of stellar mass growth is from star formation,  \n<!-- image -->  \nFIG. 11.Left panel: Ratio of the average SFR to the average baryon accretion rate ( fbdM / dt ) for halos as a function of halo mass. Right panel: same, except as a function of halo mass at z = 0 (i.e., ratio of SFR to baryon accretion rate for the progenitors of present-day halos). Lines shown the best fit model and shaded regions show the one-sigma posterior distribution.  \n<!-- image -->  \n<!-- image -->  \nFIG. 12.Left panel: Ratio of the stellar mass to halo mass at a given redshift as a function of the halo mass at z = 0. Right panel: Ratio of the stellar mass at a given redshift to the stellar mass in the descendent galaxy today, as a function of the halo mass at z = 0. Lines shown the best fit model and shaded regions show the one-sigma posterior distribution.  \n<!-- image -->  \n<!-- image -->  \nFIG. 13.Left panel: Average age of presently-existing stars in halos as a function of stellar mass. Right panel: Time required to form 50% and 90% of presently-existing stars as a function of stellar mass. Lines shown the best fit model and shaded regions show the one-sigma posterior distribution.  \n<!-- image -->  \nand hence stellar mass functions alone may be used to infer galaxy star formation rates. However, the stellar mass growth rate of high-mass ( Mh > 10 13 M /circledot ), low-redshift halos is such that it cannot be explained entirely by in situ star formation. For these halos, stellar mass functions alone cannot be used to infer their star formation rates (see also the discussion in §6 and in Appendix J). Interestingly, in the local universe the transition between merger-dominated growth and star formation growth in the local Universe is roughly the mass of the Milky Way - low mass galaxies are dominated by star formation, and galaxies in halos more massive than 10 12 M /circledot are dominated by merging.", '5.4. Baryon Conversion Efficiencies': 'Assuming that the baryon accretion rate for a halo is equal to the universal baryon fraction ( fb , which is 0.17 in our cosmology) times the overall halo mass accretion rate, we can calculate the instantaneous baryon conversion efficiency implied by our models. In the left-hand panel of Fig. 11, we show the instantaneous conversion efficiency at fixed halo mass. Remarkably, this efficiency is always highest for Milky-Way sized halos (10 12 M /circledot ), corresponding well with the peak in the stellar mass to halo mass ratio seen in Fig. 7. Moreover, this efficiency (20-40%) is constant to within a factor of two over a remarkably large redshift range, suggesting that the efficiency is only a weak function of accretion rate in halos of this mass. For more massive halos, there is a decrease in the baryon conversion efficiency following redshift 2-3. This suggests that the ability of accreted material to cool onto star-forming regions becomes impaired at lower redshifts for such halos. However, this decrease is gradual, taking place over many Gyr (see also Behroozi et al. 2013a). Rather than an abrupt change in the character of infalling gas, this may suggest that the decreasing density and accretion rate of gas make it gradually more difficult for cold clumps to form, especially in the presence of an active galactic nucleus.  \nOn the right-hand side of Fig. 11, we show the historical baryon conversion efficiency for progenitors of z = 0 halos. For massive halos at z = 0, their conversion efficiency climbs steeply towards higher redshifts as their halo mass falls to roughly 10 12 M /circledot ; once their halo mass drops below that value, they then become less and less efficient at very high redshifts. For less-massive halos, such as those that are 10 12 M /circledot at z = 0, their conversion efficiency has been increasing from early times to the present day.  \nOn the left-hand side of Fig. 12, we show extrapolated historical stellar mass to halo mass ratios-proportional to integrated baryon conversion efficiencies-similar to the bottom panel of Fig. 7. Notably, stellar mass to halo mass ratios also peak when halos reach 10 12 M /circledot , indicative of steeply falling star formation efficiencies at higher and lower masses. It would appear that the maximum integrated stellar mass efficiency is around 20-40% at all redshifts.', '5.5. Stellar Ages': "In the right-hand panel of Fig. 12, we show the historical stellar mass in progenitors of halos at z = 0 relative to the present-day stellar mass. This should not be confused with Fig. 9, which shows what amount of the currently-remaining stellar mass was in place at a given redshift; these differ mainly because of passive stellar evolution (massive stars that formed sufficiently long ago will burn out by the present day).  \nThese data allow us to derive stellar mass-weighted ages  \nas a function of stellar mass, shown in the left-hand panel of Fig. 13. These ages are consistent with the stellar populations in massive galaxies (and halos) forming at very early times and having little star formation continuing to the present day. Less-massive galaxies have younger average stellar ages, consistent with the ongoing star formation seen in the star formation histories for such galaxies in Fig. 6. While we have little information on the progenitors of galaxies less massive than 10 9 M /circledot , there is evidence that the average ages of the stellar populations may increase for such galaxies. Specifically, as shown in the right-hand panel of Fig. 13, the time required to form 50% and 90% of the galaxy's stars increases towards lower stellar masses, indicating that the star formation histories become more and more flat (constant), on average. A perfectly flat star formation history would result in an average age of about half the age of the universe, or roughly 6.7 Gyr; this is somewhat higher than that seen for 10 9 M /circledot galaxies. On the other hand, for individual dwarf galaxies, the star formation history is likely to be stochastic, which may result in a large scatter around this average value.  \nFor the most massive galaxies today ( ∼ 10 12 M /circledot ), the time required to form 50% of their stars was extremely shortonly 1-2 Gyr. This in turn translates to a high star formation rate, on the order of 200-1000 M /circledot yr -1 , which is consistent with observations of ultra-luminous infrared galaxies (ULIRGs; Magnelli et al. 2012; Michałowski et al. 2010a,b; Daddi et al. 2005; Chapman et al. 2004).", '5.6. Comparison to Other Results': "We show a comparison of our best-fit results for the stellar mass to halo mass ratio at z = 0 . 1 to previously-published results in Fig. 14, and we show comparisons at z = 1 . 0 and z = 3 . 0 in Fig. 15. Where possible, conversions to our assumed cosmology and halo mass definition have been applied. The results in this work are almost identical to our previous results in Behroozi et al. (2010), with the exception of a deviation at low halo masses due to the updated stellar mass functions used. We also compare to the new constraint of Reddick et al. (2012), which uses additional input from the correlation function and conditional stellar mass function as measured by SDSS. There is a slight discrepancy between the two results for small masses due to the different satellite fractions obtained using M peak and v peak max (still within our systematic errors), but the Reddick et al. (2012) result is only well constrained by additional data in the range 10 12 -10 14 M /circledot .  \nComparisons to other data sets at z ≤ 1 . 0 (Moster et al. 2010; Guo et al. 2010; Wang & Jing 2010; Zheng et al. 2007a; Yang et al. 2009a; Hansen et al. 2009; Lin & Mohr 2004) are discussed extensively in Behroozi et al. (2010), so we do not repeat that discussion here. The z = 3 comparison is notable because it illustrates the large discrepancies that can occur when different stellar mass functions are used at z /greaterorsimilar 3. For example, Yang et al. (2012) perform modeling for two separate highz stellar mass functions (Pérez-González et al. 2008 and Drory et al. 2005). Part of the discrepancy between the two sets of results may be due to a somewhat restrictive redshift fit, and part of it may be due to systematic biases in stellar masses at high redshifts.  \nIn Fig. 16, we show a comparison of galaxy stellar mass histories between our results and those of Leitner (2012). A direct comparison is difficult because Leitner (2012) only considers star-forming galaxies. However, according to Eq. 8, the active fraction for 10 9 . 5 M /circledot galaxies is 90%, and the ac-  \nFIG. 14.- Comparison of our best-fit model at z = 0 . 1 to previously published results. Results compared include those from our previous work (Behroozi et al. 2010), from abundance matching (Moster et al. 2013; Reddick et al. 2012; Moster et al. 2010; Guo et al. 2010; Wang & Jing 2010), from HOD/CLF modeling (Zheng et al. 2007a; Yang et al. 2012), and from cluster catalogs (Yang et al. 2009a; Hansen et al. 2009; Lin & Mohr 2004). Grey shaded regions correspond to the 68% confidence contours of Behroozi et al. (2010). The one-sigma posterior distribution for our model is shown by the red error bars.  \n<!-- image -->  \n<!-- image -->  \nFIG. 15.- Comparison of our best-fit model at z = 1 . 0 and z = 3 . 0 to previously published results. Results compared include those from our previous work (Behroozi et al. 2010), from abundance matching (Moster et al. 2013, 2010; Conroy & Wechsler 2009; Wang & Jing 2010), and from HOD/CLF modeling (Zheng et al. 2007a; Yang et al. 2012; Wake et al. 2011). Yang et al. (2012) reports best fits for two separate stellar mass functions, and we include both at z = 3 . 0. Grey shaded regions correspond to the 68% confidence contours of Behroozi et al. (2010).  \n<!-- image -->  \nFIG. 17.- Allowed range for the completeness (i.e., probability of galaxy detection) of the surveys we use as a function of redshift; the black line shows the median value, and the grey bands show 68% confidence contours. This average applies only to galaxies above authors' stated completeness limits. No strong evidence is seen for missing large numbers of high-redshift galaxies due to burstiness or dustiness.  \n<!-- image -->  \n<!-- image -->  \nz  \nFIG. 16.- Comparison of our best-fit model's stellar mass histories to those of Leitner (2012) (lines; shaded regions show the one-sigma posterior distribution of our model). Note that the Leitner (2012) results only consider star forming galaxies; ours show the average histories for all galaxies.  \ntive fraction for 10 10 M /circledot galaxies is 65%-although most of the passive 10 10 M /circledot galaxies have only recently become passive. Comparison of the stellar mass histories for those galaxies is thus somewhat more feasible. Our stellar mass histories agree remarkably well with those in Leitner (2012) except at recent times for 10 10 M /circledot galaxies (where our model includes more passive galaxies) and at early times for 10 9 . 5 M /circledot galaxies, where observations do not constrain properties of their progenitors.", '5.7. Systematic Uncertainties': 'Several aspects of the allowed parameter range are notable. Most importantly, there is no necessity for a large number of galaxies to be missed (due to dustiness or burstiness) at high redshift; indeed, the allowed incompleteness is in general less than 30%, as shown in the top panel of Fig. 17. This would imply that most high-redshift surveys are not missing a large fraction (>50%) of galaxies above their nominal detection thresholds.  \nThe allowed parameter range is consistent with having both µ = 0 and κ = 0 at z = 0; i.e., no systematic offsets in stellar mass are necessary at z = 0. This is partially by design: the parametrization of the SMHM was chosen so that it could fit the z = 0 constraints without need for µ or κ ; moreover, the cosmic and specific star formation rates are derivative constraints, meaning that they only constrain how the SMHM evolves with time. While the main fits favor some redshift evolution in the systematic parameters, it is still possible to reasonably fit all observations while disallowing all system-  \n<!-- image -->  \natic corrections (i.e., assuming that the stellar masses and star formation rates used in this paper are on average correct and that all the surveys were 100% complete). Details of this fit are presented in Appendix J. We also present discussions therein of how the allowed parameter space varies when certain data sets are not used to constrain the fit; specifically, excluding z > 6 data; excluding cosmic star formation constraints; and excluding all star formation constraints.', '5.8. Fits to Star Formation Histories': "The most commonly used star formation history fitting formula is a declining exponential with time ( SFH ( t ) ∝ e -t / τ ), used by almost all the observational papers cited in this work (see references in Tables 3 and 5). Various modifications to this form have been suggested (e.g., te -t / τ or e + t / τ in Maraston et al. 2010), and there has been increasing recent evidence that increasing star formation histories are more appropriate for massive galaxies at high redshifts.  \nOur best constraints in this study are on average star formation histories for halos at a given mass. Thus, we can place the strongest priors on the expected form of stacked galaxy star formation histories. However, we can also place some weaker constraints on what individual galaxy star formation histories should look like.  \nIndividual galaxies' histories can differ from the average star formation histories for three reasons. The first is stochasticity in the star formation rate on short timescales relative to the dynamical time of the halo, due to feedback effects as well as shot noise in the star formation rate. This variation cannot be modeled with our current approach, so we caution that the constraints on individual galaxy histories that we derive should be considered as smoothed over a period of at least 500 Myr. If enough color bands are available to avoid degeneracies with other fitting parameters, observers may consider adding the possibility of recent starbursts to their models to account for the fact that they can significantly alter observed galaxy colors.  \nThe next reason for deviations is sustained deficits in the star formation rate for satellite galaxies compared to field galaxies (e.g. Wetzel et al. 2012). As shown in Fig. 24 (see also e.g. Kravtsov et al. 2004), the satellite fraction of halos falls both with increasing mass and increasing redshift. This is not a practical concern for z > 3, where the satellite  \n<!-- image -->  \nFIG. 18.- Featureless, rising power law star formation histories are appropriate for z > 3; however, for z < 3, galaxies have peaks in their star formation rates that depend on halo mass. These figures show constraints on individual star formation histories of galaxies in halos of mass 10 11 -10 14 M /circledot as a function of time from the beginning of the universe to z = 0 . 1. These halo masses correspond to stellar masses of ∼ 10 9 , 3 × 10 10 M /circledot , 10 11 M /circledot , and 2 × 10 11 M /circledot at z = 0, respectively.  \n<!-- image -->  \n<!-- image -->  \nz  \nFIG. 19.Left panel: Median mass accretion histories for halos in narrow mass bins (0.25 dex). The shaded regions contain 68% of the spread in mass accretion histories for different halos. Right panel: Constraints on individual star formation histories for galaxies in 10 12 M /circledot halos, compared with fits from Eq. 22-24.  \n<!-- image -->  \nfraction is less than 20% for all currently-observable galaxies ( M ∗ > 10 9 . 5 M /circledot ). At lower redshifts, the observation that satellites cease to form stars relatively quickly after accretion (Wetzel et al. 2012) means that declining star formation histories may continue to be reasonable fits for satellite galaxies, even though most other galaxies of the same mass may be better fit by rising star formation histories.  \nFinally, we note that differences can arise in individual galaxy star formation histories because of the scatter in stellar masses at a given halo mass, as well as the scatter in mass accretion histories for halos (see Fig. 19). We can attempt to model these effects by sampling random mass accretion histories from the Bolshoi simulation and then sampling random stellar mass growth histories as allowed by the scatter in stellar mass at fixed halo mass (i.e., choosing random stellar mass offsets at z =8, z =1 . 0, and z =0 and using spline interpolation at intermediate times). Our results, expressed as a function of time since the beginning of the universe, are shown in Fig. 18.  \nThese results have dramatically larger error bars than those in Fig. 6 on account of the spread in mass accretion histories for halos. Nonetheless, some basic conclusions can be drawn. Star formation histories for z > 3 galaxies increase with time. Although there are different ways to parametrize these histories, a straightforward one is a direct power law form:  \nSFH ( t ) = At B (for z > 3) . (22)  \nThis is also equally capable of fitting the average star formation rates in Fig. 6 at z > 3. Our result is similar to that of Papovich et al. (2011); however, we find steeper slopes ( B ∼ 3 -4) than theirs ( B ∼ 1 . 7) because they ignore all effects of mergers that occur from z = 8 to z = 3.  \nAt lower redshifts, there is a mass-dependent turnover after which the star formation history begins to decline. This happens at z ∼ 2 -3 for 10 11 M /circledot galaxies, making a declining exponential at z = 0 a reasonable fit. However, it happens as late as z =0 . 5 for 10 9 . 5 M /circledot galaxies-meaning that a declining exponential is never a good fit for these galaxies unless they are satellites. The best fit in general for all the constraints on individual histories presented here is a double power law:  \nSFH ( t ) = A [ ( t τ ) B + ( t τ ) -C ] -1 ; (23)  \ntypical values of B and C range from 1 to 5. It may often be the case that the available data is insufficient to constrain all three shape parameters. In this case, a hybrid of the exponential decline and power law rise still provides a reasonable fit:  \nSFH ( t ) = At B exp( -t /τ ) . (24)", '6. DISCUSSION': "The existence of a 'cold mode' of gas accretion that allows efficient star formation at high redshifts and shuts off for  \nmassive galaxies past z ∼ 2 has been predicted in hydrodynamical simulations for the past decade (Birnboim & Dekel 2003). In this work, we observe a transition after z ∼ 3, where galaxies in massive halos become progressively less efficient at forming stars, even after correcting for slowing accretion rates since that time (Fig. 11, right panel). At fixed halo mass (Fig. 11, left panel), the transition is less severe, but there has been a clear reduction in the conversion efficiency since z ∼ 3 to the present day. This suggests a picture in which dense infalling gas can cool to the galaxy disk efficiently at high redshifts; however, at lower redshifts, the infalling gas becomes sparse enough that it is susceptible either to shock heating or feedback from an active galactic nucleus.  \nAnother important feature of our results is that merging galaxies' stars almost never reach the central galaxy in massive halos ( < 10% since z = 1 for halos larger than 10 13 M /circledot ). This claim is not new (Purcell et al. 2007; Conroy et al. 2007), but it is reassuring that we duplicate this conclusion with high confidence. This is in stark contrast to what is asserted in (Moster et al. 2013) (i.e., 80% of merging stars reach the central galaxy); however, their assertion is based on simulations of Milky-Way sized halos, for which we find that indeed ∼ 70% of merging galaxies' stars reach the central galaxy. Thus, the fraction of merging stars that reach the central galaxy is likely a strong function of halo mass.  \nIt is worth mentioning that the SFRs derived in Moster et al. (2013) appear reasonable despite this issue. We note that Moster et al. (2013) assumes that satellite stellar mass is fixed at the epoch of accretion, which results in less stellar mass in satellites than in this work. Additionally, Moster et al. (2013) assumes a higher value of the scatter in observed stellar mass vs. true stellar mass ( σ ) than we do, which comes from comparing technique-to-technique scatter in recovering observed stellar masses. This results in a much larger allowed evolution for massive galaxies between z = 1 and z = 0, which is very sensitive to the value of σ (Behroozi et al. 2010). The combination of these two assumptions (less stellar mass in satellites and greater allowed evolution in the stellar mass of massive galaxies) roughly cancels the effect on the SFR of assuming that more stellar mass in satellites is allowed to merge into the central galaxy. We note, however, that our satellite stellar masses are already on the low side of what is allowed by clustering constraints (Reddick et al. 2012), leaving little room to reduce them further as is done in the Moster et al. (2013) analysis. In addition, the scatter in the difference between two stellar mass estimation techniques tends to overestimate the scatter between any one technique and the true value by a factor of ∼ √ 2. This factor is consistent with other independent estimates of σ , for example, by evaluating the posterior distribution of allowed recovered stellar masses (Kajisawa et al. 2009; Conroy et al. 2009). Thus, both the actual evolution in the stellar mass of massive galaxies and the fraction of stars in mergers which reach the central galaxy may not be as large the values implied by the Moster et al. (2013) analysis.  \nWe note further that the amount of stellar mass that was contained in disrupted galaxies can easily reach five times the amount of stellar mass in the central galaxy for high-mass halos (Fig. 9), consistent with observations in Gonzalez et al. (2005). Regardless of uncertainties in determining stellar masses, there is no possible way to funnel any significant fraction of this mass into the central galaxy without badly violating observed limits on the number density of massive galaxies. Rather, stars from merging galaxies must end up in a  \nmore extended distribution to avoid confusion with the central galaxy, a topic that will be explored in depth in a future paper (Behroozi et al., in prep). From the opposite perspective, galaxy stellar populations contain mostly stars from the most-massive progenitor, with some contamination from major mergers. Surveys of stars within a galaxy's optical radius (such as for the Milky Way) thus capture a potentially limited range of different star-forming environments, at least for stars formed since z = 1. A wider range of environments can be sampled by considering the ICL and halo stars in general.  \nThe most efficient halo mass for converting baryons into stars appears to be close to 10 12 M /circledot , regardless of redshift. The progenitors of today's massive halos ( /greaterorsimilar 10 14 M /circledot ) passed through this mass quickly, and as a result formed most of their stars in a short time: a mere 1-2 Gyr (Fig. 13). Given that baryon to star conversion becomes so inefficient above and below this halo mass, much of the shape of the star formation histories seen in this paper (Fig. 6) could be explained by considering the duration of time that a given halo's mass accretion history spends close to this most efficient halo mass. That is to say, the peakiness of the star formation history depends on how quickly the halo mass accretion history passes through halos of mass ∼ 10 12 M /circledot . As lower-mass halos have slower relative accretion rates, they have flatter star formation histories; e.g., Milky Way-sized halos have lingered in this mass range from z ∼ 1, and so have the current highest stellar mass to halo mass ratio, as well as star formation histories that track their halo mass accretion rates well at late times.  \nLow-mass halos (10 11 M /circledot ) show an unusual behavior in terms of their baryon conversion efficiency. In the left panel of Fig. 11, it would appear that the instantaneous baryon conversion efficiency of 10 11 M /circledot halos is high at high redshifts ( z > 5), low at intermediate redshifts (2 < z < 5), and then high again at low redshifts. While this behavior matches what is seen in the integrated baryon conversion efficiency for such halos (Fig. 12, left panel), the corresponding galaxies are often at the faint edge of the completeness limit for stellar mass surveys. For that reason, observational biases may be the most likely source of this behavior.  \nOur results have importance for future observational studies. The current levels of systematic errors in the inference of stellar masses and star formation rates are the dominant systematic uncertainty in this work. As discussed in §3, stellar mass functions at z > 3 can disagree by up to 0.5 dex, well beyond any errors expected from sample variance. Star formation rates are equally concerning, both at high redshifts and at low redshifts, especially in terms of specific star formation rates. As a result, systematic errors are the single most important aspect preventing better understanding of the formation histories of massive galaxies. As discussed in §5.8, we have provided new fitting formulas for galaxy star formation histories, which should help reduce one of the largest sources of systematic error especially for low-mass and high-redshift galaxies. New studies of stellar mass functions and star formation rates from deep multiwavelength surveys like CANDELS (Koekemoer et al. 2011; Grogin et al. 2011), analyzed with self-consistent star formation histories and assumptions about systematic errors, will provide a large step forward in constraining the model presented here.", '7. CONCLUSIONS': "Wehavepresented a comprehensive analysis of average star formation rates and histories in galaxies from z = 0 to z = 8,  \nand their connection to the underlying growth and merging of dark matter halos, along with a treatment of the inherent uncertainties. Our approach provides a self-consisent picture of the growth of galaxies over this epoch in the cosmological context of Λ CDM structure formation. The model is able to match the evolution of cosmic star formation rates, specific star formation rates, and the stellar mass function of galaxies over the last 13.2 billion years of cosmic time. Our main findings are as follows:  \n- 1. Halos of mass ∼ 10 12 M /circledot appear to be the most efficient at forming stars at every epoch, with baryon conversion efficiencies of 20-40% over nearly the entire redshift range of this study. Halos at higher and lower masses are less efficient by orders of magnitude, especially at low redshifts.\n- 2. The baryon conversion efficiency of more massive halos is still reasonably high (10-20%) until z =2 -3, after which time it takes a steep downturn. This is consistent with expectations of reduced cold mode accretion onto the corresponding galaxies.\n- 3. We have characterized the fraction of stars in galaxies due to merging vs star formation. At z > 1, galaxy buildup in all but the most massive galaxies was strongly dominated by star formation. At the present day, the transition between merger-dominated buildup (at high mass) and star formation-dominated buildup (at low mass) occurs at roughly the mass of the Milky Way.\n- 4. We confirm previous results in support of most merging galaxies' stars being disrupted in the ICL for massive halos. We find, however, that the fraction of merging stars that reach the central galaxy may increase strongly with decreasing halo mass, to about 70% for Milky Way-sized halos.\n- 5. For massive galaxies ( M ∗ > 10 11 M /circledot ), we predict that at least as much stellar mass is present in the ICL as in the galaxy itself out to intermediate redshifts, possibly as high as z ∼ 1 for the largest clusters (although this depends on the exact cut between ICL and galaxy stellar mass).\n- 6. Star formation histories peak at z ∼ 3 for massive galaxies and are increasing nearly to the present day for small galaxies. We have presented new fitting formulas for galaxy star formation histories: a featureless power law for high redshifts (Eq. 22), and more general formulas valid for all redshifts (Eq. 23 and 24).\n- 7. Upturns in the stellar mass function at low stellar masses mean that dwarf galaxies ( M ∗ < 10 8 . 5 M /circledot ) have\n- a higher stellar mass to halo mass ratio than that predicted by extrapolating results from more massive galaxies. The upturn is most robustly observed at low redshifts, but there are also indications at higher redshifts for a steepening slope to the SMF at low stellar masses. The stellar mass to halo mass relation thus cannot be fit with a double power law alone.\n- 8. Systematic errors remain at high levels ( ∼ 0 . 2 -0 . 3 dex in stellar masses, 0 . 15 -0 . 3 dex in SFRs), and disagreement between galaxy abundances for high-redshift surveys ( z > 3) is substantial.\n- 9. Galaxies with stellar masses below 10 10 M /circledot have the least well-constrained properties in this study. This could be improved by future surveys and by incorporating information on the star formation histories of local galaxies.  \nSupport for this work was provided by an HST Theory grant. Program number HST-AR-12159.01-A was provided by NASA through a grant from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Incorporated, under NASA contract NAS5-26555. This research was also supported in part by the National Science Foundation under Grant No. NSF PHY11-25915, through a grant to KITP during the program 'First Galaxies and Faint Dwarfs', and through a ARCS Foundation Graduate Fellowship to PSB.  \nWethank Anatoly Klypin and Joel Primack for access to the Bolshoi and MultiDark simulations and Michael Busha and the LasDamas collaboration for access to the Consuelo simulation. The latter was run using computational resources at SLAC. We thank John Moustakas for access to the PRIMUS stellar mass functions prior to publication. This work has benefitted from many interesting conversations with Tom Abel, Steve Allen, Andrew Benson, Andreas Berlind, Kevin Bundy, Joanne Cohn, Sandy Faber, Dusan Keres, David Koo, Andrey Kravtsov, Alexie Leauthaud, Avi Loeb, Yu Lu, ChungPei Ma, John Moustakas, Joel Primack, Rachel Reddick, Beth Reid, Daniel Stark, Chuck Steidel, Tommaso Treu, David Wake, Andrew Wetzel, Frank van den Bosch, Alexey Vikhlinin, and Naoki Yoshida, We thank Yu Lu, Ari Maller, Simone Wienmann, and Anatoly Klypin for helpful comments on an earlier draft, and Surhud More and Michal Michałowski for catching typographical mistakes. We additionally the anonymous referee for many suggestions which improved this paper.", 'REFERENCES': "Berrier, J. C., Bullock, J. S., Barton, E. J., Guenther, H. D., Zentner, A. R., &Wechsler, R. H. 2006, ApJ, 652, 56 Béthermin, M., Doré, O., & Lagache, G. 2012, A&A, 537, L5 Birnboim, Y., & Dekel, A. 2003, MNRAS, 345, 349 Blanton, M. R., & Roweis, S. 2007, AJ, 133, 734 Borgani, S., & Kravtsov, A. 2011, Advanced Science Letters, 4, 204 Bouwens, R. J., et al. 2012, ApJ, 754, 83 -. 2011, ApJ, 737, 90  \nBoylan-Kolchin, M., Bullock, J. S., & Kaplinghat, M. 2012, MNRAS, 422, 1203  \nBaldry, I. K., Glazebrook, K., & Driver, S. P. 2008, MNRAS, 388, 945  \nBehroozi, P. S., Conroy, C., & Wechsler, R. H. 2010, ApJ, 717, 379  \nBehroozi, P. S., Loeb, A., & Wechsler, R. H. 2012, arXiv:1208.0334  \nBehroozi, P. S., Wechsler, R. H., & Conroy, C. 2013a, ApJ, 762, L31 Behroozi, P. S., Wechsler, R. H., & Wu, H.-Y. 2013b, ApJ, 762, 109 Behroozi, P. S., Wechsler, R. H., Wu, H.-Y., Busha, M. T., Klypin, A. A., & Primack, J. R. 2013c, ApJ, 763, 18 Berlind, A. A., & Weinberg, D. H. 2002, ApJ, 575, 587 Bernardi, M., Shankar, F., Hyde, J. B., Mei, S., Marulli, F., & Sheth, R. K.  \n2010, MNRAS, 404, 2087  \n- Bradley, L. D., et al. 2012, ApJ, 760, 108  \nBrammer, G. B., et al. 2011, ApJ, 739, 24  \n- Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000\n- Bryan, G. L., & Norman, M. L. 1998, ApJ, 495, 80\n- Bullock, J. S., Wechsler, R. H., & Somerville, R. S. 2002, MNRAS, 329, 246  \nCaputi, K. I., Cirasuolo, M., Dunlop, J. S., McLure, R. J., Farrah, D., &  \nAlmaini, O. 2011, MNRAS, 413, 162  \n- Chabrier, G. 2003, Publications of the Astronomical Society of the Pacific, 115, 763\n- Chapman, S. C., Smail, I., Windhorst, R., Muxlow, T., & Ivison, R. J. 2004, ApJ, 611, 732\n- Colín, P., Klypin, A. A., Kravtsov, A. V., & Khokhlov, A. M. 1999, ApJ, 523, 32\n- Conroy, C., & Gunn, J. E. 2010, ApJ, 712, 833\n- Conroy, C., Gunn, J. E., & White, M. 2009, ApJ, 699, 486\n- Conroy, C., & Wechsler, R. H. 2009, ApJ, 696, 620\n- Conroy, C., Wechsler, R. H., & Kravtsov, A. V. 2006, ApJ, 647, 201\n- Conroy, C., Wechsler, R. H., & Kravtsov, A. V. 2007, ApJ, 668, 826\n- Conroy, C., White, M., & Gunn, J. E. 2010, ApJ, 708, 58\n- Crocce, M., Pueblas, S., & Scoccimarro, R. 2006, MNRAS, 373, 369\n- Cucciati, O., et al. 2012, A&A, 539, A31\n- Daddi, E., et al. 2005, ApJ, 631, L13  \n-. 2007, ApJ, 670, 156  \n- Donahue, M., et al. 2010, ApJ, 715, 881\n- Drory, N., et al. 2009, ApJ, 707, 1595\n- Drory, N., Salvato, M., Gabasch, A., Bender, R., Hopp, U., Feulner, G., & Pannella, M. 2005, ApJ, 619, L131\n- Dunkley, J., Bucher, M., Ferreira, P. G., Moodley, K., & Skordis, C. 2005, MNRAS, 356, 925\n- Dunne, L., et al. 2009, MNRAS, 394, 3\n- Eddington, A. S. 1913, MNRAS, 73, 359\n- Eddington, Sir, A. S. 1940, MNRAS, 100, 354\n- Fakhouri, O., Ma, C.-P., & Boylan-Kolchin, M. 2010, MNRAS, 406, 2267 Feulner, G., Gabasch, A., Goranova, Y., Hopp, U., & Bender, R. 2008, in Relativistic Astrophysics Legacy and Cosmology - Einstein's, ed. B. Aschenbach, V. Burwitz, G. Hasinger, & B. Leibundgut, 310-+  \nFirmani, C., & Avila-Reese, V. 2010, ApJ, 723, 755  \n- Gavazzi, R., Treu, T., Rhodes, J. D., Koopmans, L. V. E., Bolton, A. S., Burles, S., Massey, R. J., & Moustakas, L. A. 2007, ApJ, 667, 176 Gonzalez, A. H., Zabludoff, A. I., & Zaritsky, D. 2005, ApJ, 618, 195\n- Gonzalez, V., Bouwens, R., llingworth, G., Labbe, I., Oesch, P., Franx, M., & Magee, D. 2012, arXiv:1208.4362\n- González, V., Labbé, I., Bouwens, R. J., Illingworth, G., Franx, M., & Kriek, M. 2011, ApJ, 735, L34\n- Grazian, A., et al. 2011, A&A, 532, A33  \nGrogin, N. A., et al. 2011, ApJS, 197, 35  \n- Guo, Q., White, S., Li, C., & Boylan-Kolchin, M. 2010, MNRAS, 404, 1111 Haario, H., Saksman, E., & Tamminen, J. 2001, Bernoulli, 7, pp. 223\n- Hammer, F., Puech, M., Chemin, L., Flores, H., & Lehnert, M. D. 2007, ApJ, 662, 322\n- Hansen, S. M., Sheldon, E. S., Wechsler, R. H., & Koester, B. P. 2009, ApJ, 699, 1333\n- Hopkins, A. M., & Beacom, J. F. 2006, ApJ, 651, 142  \nJenkins, A. 2010, MNRAS, 403, 1859  \n- Kajisawa, M., et al. 2009, ApJ, 702, 1393  \nKajisawa, M., Ichikawa, T., Yamada, T., Uchimoto, Y. K., Yoshikawa, T.,  \nAkiyama, M., & Onodera, M. 2010, ApJ, 723, 129  \nKarim, A., et al. 2011, ApJ, 730, 61  \n- Kistler, M. D., Yüksel, H., Beacom, J. F., Hopkins, A. M., & Wyithe, J. S. B. 2009, ApJ, 705, L104  \nKlypin, A. A., Trujillo-Gomez, S., & Primack, J. 2011, ApJ, 740, 102  \n- Koekemoer, A. M., et al. 2011, ApJS, 197, 36\n- Komatsu, E., et al. 2009, ApJS, 180, 330  \n-. 2011, ApJS, 192, 18  \n- Kravtsov, A., & Klypin, A. 1999, ApJ, 520, 437\n- Kravtsov, A. V., Berlind, A. A., Wechsler, R. H., Klypin, A. A., Gottloeber, S., Allgood, B., & Primack, J. R. 2004, ApJ, 609, 35  \nKravtsov, A. V., Klypin, A. A., & Khokhlov, A. M. 1997, ApJ, 111, 73  \n- Labbé, I., et al. 2010, ApJ, 716, L103\n- Labbe, I., et al. 2012, arXiv:1209.3037\n- Le Borgne, D., Elbaz, D., Ocvirk, P., & Pichon, C. 2009, A&A, 504, 727 Leauthaud, A., et al. 2012, ApJ, 744, 159  \nLee, K.-S., et al. 2011, ApJ, 733, 99  \n-. 2012, ApJ, 752, 66  \n- Lee, K.-S., Giavalisco, M., Conroy, C., Wechsler, R. H., Ferguson, H. C., Somerville, R. S., Dickinson, M. E., & Urry, C. M. 2009, ApJ, 695, 368 Leitner, S. N. 2012, ApJ, 745, 149  \n- Lin, Y.-T., & Mohr, J. J. 2004, ApJ, 617, 879\n- Lu, Y., Mo, H. J., Katz, N., & Weinberg, M. D. 2012, MNRAS, 421, 1779 Lucatello, S., Gratton, R. G., Beers, T. C., & Carretta, E. 2005, ApJ, 625, 833\n- Ly, C., Lee, J. C., Dale, D. A., Momcheva, I., Salim, S., Staudaher, S., Moore, C. A., & Finn, R. 2011a, ApJ, 726, 109  \nLy, C., Malkan, M. A., Hayashi, M., Motohara, K., Kashikawa, N.,  \nShimasaku, K., Nagao, T., & Grady, C. 2011b, ApJ, 735, 91  \nMagnelli, B., Elbaz, D., Chary, R. R., Dickinson, M., Le Borgne, D., Frayer,  \nD. T., & Willmer, C. N. A. 2011, A&A, 528, A35+  \nMagnelli, B., et al. 2012, A&A, 539, A155  \n- Mandelbaum, R., Seljak, U., Kauffmann, G., Hirata, C. M., & Brinkmann, J. 2006, MNRAS, 368, 715\n- Maraston, C., Pforr, J., Renzini, A., Daddi, E., Dickinson, M., Cimatti, A., & Tonini, C. 2010, MNRAS, 407, 830\n- Marchesini, D., van Dokkum, P. G., Förster Schreiber, N. M., Franx, M., Labbé, I., & Wuyts, S. 2009, ApJ, 701, 1765  \nMarchesini, D., et al. 2010, ApJ, 725, 1277  \n- Marín, F. A., Wechsler, R. H., Frieman, J. A., & Nichol, R. C. 2008, ApJ, 672, 849\n- McBride, J., Fakhouri, O., & Ma, C.-P. 2009, MNRAS, 398, 1858 McLure, R. J., et al. 2011, MNRAS, 418, 2074\n- Michałowski, M., Hjorth, J., & Watson, D. 2010a, A&A, 514, A67  \nMichałowski, M. J., Watson, D., & Hjorth, J. 2010b, ApJ, 712, 942  \n- More, S., van den Bosch, F. C., Cacciato, M., Mo, H. J., Yang, X., & Li, R. 2009, MNRAS, 392, 801\n- Mortlock, A., Conselice, C. J., Bluck, A. F. L., Bauer, A. E., Grützbauch, R., Buitrago, F., & Ownsworth, J. 2011, MNRAS, 413, 2845\n- Moster, B. P., Naab, T., & White, S. D. M. 2013, MNRAS, 428, 3121  \nMoster, B. P., Somerville, R. S., Maulbetsch, C., van den Bosch, F. C.,  \nMacciò, A. V., Naab, T., & Oser, L. 2010, ApJ, 710, 903  \nMoustakas, J., et al. 2013, arXiv:1301.1688  \n- Muzzin, A., Marchesini, D., van Dokkum, P. G., Labbé, I., Kriek, M., & Franx, M. 2009, ApJ, 701, 1839\n- Nagamine, K., Ostriker, J. P., Fukugita, M., & Cen, R. 2006, ApJ, 653, 881 Neyrinck, M. C., Hamilton, A. J. S., & Gnedin, N. Y. 2004, MNRAS, 348, 1 Noeske, K. G., et al. 2007, ApJ, 660, L43\n- Papovich, C., Finkelstein, S. L., Ferguson, H. C., Lotz, J. M., & Giavalisco, M. 2011, MNRAS, 412, 1123\n- Pérez-González, P. G., et al. 2008, ApJ, 675, 234\n- Purcell, C. W., Bullock, J. S., & Zentner, A. R. 2007, ApJ, 666, 20\n- Rawle, T. D., et al. 2012, ApJ, 747, 29\n- Reddick, R. M., Wechsler, R. H., Tinker, J. L., & Behroozi, P. S. 2012, arXiv:1207.2160\n- Reddy, N. A., Pettini, M., Steidel, C. C., Shapley, A. E., Erb, D. K., & Law, D. R. 2012, ApJ, 754, 25  \nReddy, N. A., & Steidel, C. C. 2009, ApJ, 692, 778  \nRiebe, K., et al. 2011, arXiv:1109.0003  \n- Robotham, A. S. G., & Driver, S. P. 2011, MNRAS, 413, 2570\n- Rujopakarn, W., et al. 2010, ApJ, 718, 1171\n- Salim, S., et al. 2007, ApJS, 173, 267  \nSalmi, F., Daddi, E., Elbaz, D., Sargent, M. T., Dickinson, M., Renzini, A.,  \nBethermin, M., & Le Borgne, D. 2012, ApJ, 754, L14  \nSchaerer, D., & de Barros, S. 2010, A&A, 515, A73+  \n- Shankar, F., Lapi, A., Salucci, P., De Zotti, G., & Danese, L. 2006, ApJ, 643, 14\n- Shim, H., Colbert, J., Teplitz, H., Henry, A., Malkan, M., McCarthy, P., & Yan, L. 2009, ApJ, 696, 785  \nSmolˇci'c, V., et al. 2009, ApJ, 690, 610  \nSobral, D., Smail, I., Best, P. N., Geach, J. E., Matsuda, Y., Stott, J. P.,  \nCirasuolo, M., & Kurk, J. 2013, MNRAS, 428, 1128  \nSpringel, V. 2005, MNRAS, 364, 1105  \n- Stark, D. P., Ellis, R. S., Bunker, A., Bundy, K., Targett, T., Benson, A., & Lacy, M. 2009, ApJ, 697, 1493\n- Stark, D. P., Schenker, M. A., Ellis, R., Robertson, B., McLure, R., & Dunlop, J. 2013, ApJ, 763, 129\n- Tadaki, K.-I., Kodama, T., Koyama, Y., Hayashi, M., Tanaka, I., & Tokoku, C. 2011, PASJ, 63, 437\n- Tasitsiomi, A., Kravtsov, A. V., Wechsler, R. H., & Primack, J. R. 2004, ApJ, 614, 533\n- Tinker, J., Kravtsov, A. V., Klypin, A., Abazajian, K., Warren, M., Yepes, G., Gottlöber, S., & Holz, D. E. 2008, ApJ, 688, 709  \nTrenti, M., & Stiavelli, M. 2008, ApJ, 676, 767  \nTumlinson, J. 2007a, ApJ, 665, 1361  \n-. 2007b, ApJ, 664, L63  \n- Twite, J. W., Conselice, C. J., Buitrago, F., Noeske, K., Weiner, B. J., Acosta-Pulido, J. A., & Bauer, A. E. 2012, MNRAS, 420, 1061  \n- Vale, A., & Ostriker, J. P. 2004, MNRAS, 353, 189 -. 2006, MNRAS, 371, 1173 van der Burg, R. F. J., Hildebrandt, H., & Erben, T. 2010, A&A, 523, A74+ van Dokkum, P. G. 2008, ApJ, 674, 29 van Dokkum, P. G., & Conroy, C. 2010, Nature, 468, 940 Wake, D. A., et al. 2011, ApJ, 728, 46 Wang, L., et al. 2013, MNRAS Wang, L., & Jing, Y. P. 2010, MNRAS, 402, 1796 Watson, D. F., Berlind, A. A., & Zentner, A. R. 2012, ApJ, 754, 90 Wechsler, R. H., Bullock, J. S., Primack, J. R., Kravtsov, A. V., & Dekel, A. 2002, ApJ, 568, 52 Weinmann, S. M., Neistein, E., & Dekel, A. 2011, MNRAS, 417, 2737 Weinmann, S. M., Pasquali, A., Oppenheimer, B. D., Finlator, K., Mendel, J. T., Crain, R. A., & Macciò, A. V. 2012, MNRAS, 426, 2797 Wetzel, A. R., Tinker, J. L., & Conroy, C. 2012, MNRAS, 424, 232 Wetzel, A. R., & White, M. 2010, MNRAS, 403, 1072 Whitaker, K. E., van Dokkum, P. G., Brammer, G., & Franx, M. 2012, ApJ, 754, L29\n- White, M., Zheng, Z., Brown, M. J. I., Dey, A., & Jannuzi, B. T. 2007, ApJ,\n- 655, L69 Wilkins, S. M., Trentham, N., & Hopkins, A. M. 2008, MNRAS, 385, 687 Yang, X., Mo, H. J., & van den Bosch, F. C. 2009a, ApJ, 695, 900 -. 2009b, ApJ, 693, 830 Yang, X., Mo, H. J., van den Bosch, F. C., Zhang, Y., & Han, J. 2012, ApJ, 752, 41 Yang, X., et al. 2003, MNRAS, 339, 1057 Yoshida, M., et al. 2006, ApJ, 653, 988 Zheng, W., et al. 2012, Nature, 489, 406 Zheng, X. Z., Bell, E. F., Papovich, C., Wolf, C., Meisenheimer, K., Rix, H.-W., Rieke, G. H., & Somerville, R. 2007a, ApJ, 661, L41  \nZheng, Z., Coil, A. L., & Zehavi, I. 2007b, ApJ, 667, 760", 'A. IMPACT OF ASSUMING AVERAGE STAR FORMATION HISTORIES ON DERIVED STAR FORMATION RATES': 'As mentioned in the introduction to §2, there is a potential error in our method because we propagate average star formation histories for halos at a given mass and redshift, as opposed to individual star formation histories. We discuss effects from both a distribution in star formation rates for galaxies at a given halo mass and from a distribution of halo mass accretion histories.  \nBimodality has been observed in the star formation rates of galaxies at fixed stellar mass, partially but not exclusively driven by quenched satellite galaxies (see, e.g., Wetzel et al. 2012). This may appear to be a problem for our method because the rate of specific mass loss (see §2.3) will be different depending on whether a given galaxy is active or passive. Remarkably, however, our method is robust to any distribution of star formation rates in galaxies of the same stellar mass. Because our method only depends on the total stellar mass buildup in halos of a given mass Mh at any given time, the fact that some galaxies get a larger share of that buildup than others is not a concern. If certain galaxies sustain a higher growth rate for a long time, then the scatter in stellar mass at fixed halo mass will certainly increase; however, the possibility of scatter changing with time is an effect that we include explicitly in our model. Bimodality in star formation rates is only a problem if it results in a distribution of stellar masses at fixed halo mass that cannot be reasonably modeled by a log-normal distribution (i.e., the distribution we assume in this work). So far, no evidence has emerged for any other distribution to be preferred (see comments in §2.2.1 where ξ is defined).  \nThere is more possibility for error because we do not account for correlations between the star formation history (SFH) and the mass accretion history of a halo. For example, a rapidly accreting halo might have a more recent star formation history than other halos at the same mass. In this case, when it transitions from one mass bin to another, we may use the incorrect SFH in Eq. B2. This effect is balanced by halos that accrete more slowly than other halos of the same mass, which cancels the effect to first order.  \nIt is worth working through a more quantitative example to verify this assertion. Consider a population of galaxies all with the same stellar mass M ∗ = 10 10 M /circledot today. Assume that each galaxy has had a constant SSFR for its entire lifetime (i.e., an exponential growth history), and moreover that the distribution of the SSFRs is log-normal with width 0.3 dex (similar to that seen in Salim et al. 2007) and median 10 -9 yr -1 . The average (as opposed to median) star formation rate for this population is 12.7 M /circledot yr -1 ; given the star formation histories just specified, this translates into an average rate of change of stellar mass of 9.3 M /circledot yr -1 over the next 100 Myr. If, on the other hand, we had simply assumed that all of the galaxies had star formation histories corresponding to a constant SSFR of 10 -9 yr -1 , we would interpret a buildup in stellar mass of 9.3 M /circledot yr -1 as implying an instantaneous star formation rate of 13.4 M /circledot yr -1 , an error of 6%. Thus, we conclude that the errors on the average quantities considered here introduced by lumping together star formation histories for halos with different mass accretion histories are strongly subdominant to existing observational uncertainties.', 'B. MERGING GALAXIESAND MASS ACCRETION': "Because of passive stellar evolution, calculating the amount of stellar mass remaining in a galaxy (required for determining the star formation rate in Eq. 1), one must know the entire star formation history. For a large simulation, keeping track of the individual stellar mass histories (covering hundreds of timesteps) for millions of galaxies repeated over several million realizations of the stellar mass - halo mass relation is inefficient. Instead, since we are interested in only the average star formation rate as a function of mass, we take an elegant (although approximate) shortcut.  \nInstead of keeping track of every halo/galaxy individually, we instead keep track of average star formation histories at each timestep in bins of halo mass (which is a rank-3 tensor, SFH t f orm t , m ; t is the time now, t f orm is the timestep at which the stellar mass was formed, and m is the halo mass). For central halos, we use the virial mass (as discussed in §4); for satellite halos, we use the mass of the largest progenitor (the peak mass) in the satellite's mass accretion history. At every timestep, we calculate the number of halos that transition between one mass bin and another ( T mp t , mn ), i.e., the number of halos that had a mass mp at timestep t and a mass mn at timestep t + 1 (see Appendices G and H). Similarly, we calculate the merger rates for subhalos at each timestep as a function of the subhalo peak mass and the target host halo mass ( M mp t , mn ), as discussed in Appendix I. Keeping track of star formation histories then becomes an exercise in tensor multiplication, which is fast regardless of the number of halos in the simulation.  \nTo convert Eq. 1 into a full tensor equation, we need a small amount of additional notation. The stellar mass at a particular timestep t and halo mass m becomes SMt , m , the star formation rate at a particular timestep and halo mass becomes SFRt , m , and the number counts at a given timestep for halos in a given mass bin becomes Nt , m . In addition, we denote the difference in stellar mass from the previous timestep by ∆ SMt , m . Thus, Eq. 1 formally becomes:  \n∆ SMt , m =   SMt , m -∑ t f orm SFH t f orm t , m (1 -C loss( t , t f orm -∆ t , t f orm ))   SFRt , m = fSFR ( m , t ) ∆ SMt , m ∆ t (1 -C loss( t , t -∆ t , t )) C loss( tnow , ta , tb ) = | tb -ta | -1 ∫ tb ta f loss( tnow -x ) dx , (B1)  \nwhere C loss( tnow , ta , tb ) is the average stellar mass loss fraction (a.k.a., 'recycling' fraction) at the present time for stars that formed at a constant rate from ta to tb . (Note that the indices t and m are not meant to be summed over).  \nTo calculate star formation histories, we need to know both internal star formation and the amount of stars accreted in mergers. As discussed in §2.4, purely definitional changes in what constitutes the boundary between a galaxy and the ICL can effectively shift stellar mass from the ICL into the galaxy. Not knowing the relative fraction of each, we assume that merging stellar mass from satellites takes priority over stellar mass from the ICL. For convenience, we define SM ' t , m and ICL ' t , m be the total stellar mass in merging satellites and the total stellar mass in the main progenitors' ICL, respectively; we also define ICLH t f orm t , m to be the star formation history of the ICL (analogously to SFH t f orm t , m ). The new star formation histories will then be the sum of the old star formation histories from all progenitors, with the addition of new stars formed at the latest timestep and material from mergers and the ICL:  \nSFH t f orm t + ∆ t , m = ∑ mp [ SFH t f orm t , mp ( T mp t , m + C merge t , m M mp t , m ) + ICLH t f orm t , mp T mp t , m C def t , m ] ( Nt + ∆ t , m ) -1 + SFRt , m δ t f orm t ∆ t (B2) C merge t , m =min ( 1 , (1 -fSFR ( m , t )) ∆ SMt + ∆ t , m SM ' t , m ) (B3) C def t , m =max ( 0 , (1 -fSFR ( m , t )) ∆ SMt + ∆ t , m -SM ' t , m ICL ' t , m ) (B4) SM ' t , m = ∑ mp SFH t f orm t , mp M mp t , m (1 -C loss( t + ∆ t , t f orm -∆ t , t f orm ))( Nt + ∆ t , m ) -1 (B5) (B6)  \nICL ' t , m = ∑ mp ICLH t f orm t , mp T mp t , m (1 -C loss( t + ∆ t , t f orm -∆ t , t f orm ))( Nt + ∆ t , m ) -1 .  \n(Again, the indices t and m are not meant to be summed over; δ t f orm t is the standard identity tensor: δ t f orm t = 1 if t = t f orm and 0 otherwise). Note that, in this present analysis, we do not separate the star formation histories of satellites and centrals. As discussed in §5.8, satellites at high redshifts ( z > 1) tend to have short infall timescales, so they have roughly similar star formation histories as centrals. At later times, this approximation is less valid; however, we can nonetheless correctly recover the average star formation history across all halos (centrals and satellites) at a given halo mass, as discussed in Appendix A.  \nCalculating the ICL formation history as well as the total amount of stars in the ICL is similar:  \nICLH t f orm t + ∆ t , m = ∑ mp ( SFH t f orm t , mp (1 -C merge t , m ) M mp t , m + ICLH t f orm t , mp T mp t , m (1 -C def t , m ) ) ( Nt + ∆ t , m ) -1 (B7) ICLt , m = ∑ t f orm ICLH t f orm t , mp (1 -C loss( t , t f orm -∆ t , t f orm )) . (B8)  \nNote that the merging satellite ICL is not deposited into the central galaxy's ICL, which is the same approach taken by Conroy et al. (2007) and Purcell et al. (2007). Because the ICL of a merging satellite halo will get stripped long before the stars in its associated galaxy, the stripped satellite's ICL will end up in a much more extended, faint distribution than the host's ICL.  \nThis formal specification is equivalent to a discrete differential equation for the stellar mass history. It is thus important to properly specify the boundary conditions; that is, the stellar mass history at t = 0. In the absence of constraints on stellar mass functions above z ≈ 9, we extrapolate the fit to the stellar mass - halo mass relation to z = 15 and assume a uniform star formation rate prior to that redshift. While incorrect in detail at z = 15, this choice has little effect on the star formation rates or histories of galaxies at later (and more visible) times.", 'C. CALCULATINGOBSERVABLES': 'Given the large number of equations in our framework, it is appropriate to describe exactly how a point in the parameter space of the cosmic stellar mass history (Eqs. 3-5) and uncertainties (Eqs. 9-18) gets converted into quantities comparable to published observations. The three observables in our case are stellar mass functions, specific star formation rates, and the cosmic star formation rate density. We deal with each in turn.', 'C.1. The Observed Stellar Mass Function': 'Equation 3 gives the median stellar mass as a function of halo mass for the entire redshift range we consider. If no scatter between stellar mass and halo mass were present, the conversion from the halo mass function to the fiducial stellar mass function ( φ f id ) would be straightforward application of the chain rule:  \nφ fid( M ∗ ) = dN d log 10 M ∗ = dN d log 10 Mh ( M ∗ ) ( d log 10 M ∗ d log 10 Mh ) -1 . (C1)  \nIn other words, the stellar mass function is equal to the halo mass function divided by the logarithmic derivative of Eq. 3. As discussed in Behroozi et al. (2010), this may be convolved with the scatter in observed stellar mass at fixed halo mass ( P ( ∆ M ∗ )) to give the true stellar mass function:  \nφ true( M ∗ ) = ∫ ∞ -M ∗ φ fid( M ∗ + ∆ M ∗ ) P ( ∆ M ∗ ) d ∆ M ∗ . (C2)  \nIn our model, both the scatter in true stellar mass at fixed halo mass (Eq. 5) and the scatter in observed stellar mass at fixed true stellar mass (Eq. 11) are log-normal distributions; hence P ( ∆ M ∗ ) is a log-normal distribution with width equal to the sum in quadrature of the widths of the two sources of scatter ( σ 2 + ξ 2 ).  \n√ To convert from the true stellar mass function to the observed stellar mass function, we need to correct not only for observational completeness (Eq. 13), but also for systematic biases in the measured stellar masses (Eqs. 9 and 10). This latter correction is made somewhat trickier by the fact that we only know the passive fraction as a function of the measured stellar mass. The true passive fraction, f passive , true( M ∗ ) is related to the observed passive fraction by:  \nf passive , obs( M ∗ , meas) = φ true( M ∗ , meas10 -µ ) f passive , true( M ∗ , meas10 -µ ) φ true( M ∗ , meas10 -µ ) f passive , true( M ∗ , meas10 -µ ) + φ true( M ∗ , meas10 -µ -κ )(1 -f passive , true( M ∗ , meas10 -µ -κ )) . (C3)  \nIn the case that κ =0, then this equation reduces to f passive , true( M ∗ , meas10 -µ ) = f passive , obs( M ∗ , meas). Otherwise, this equation becomes a recurrence relation:  \nf passive , true( M ∗ , meas10 -µ ) 1 -f passive , true( M ∗ , meas10 -µ -κ ) = φ true( M ∗ , meas10 -µ -κ ) φ true( M ∗ , meas10 -µ ) [ (1 -f passive , obs( M ∗ , meas)) -1 -1 ] . (C4)  \nWe therefore compromise by decoupling the true passive fraction from the value of κ . We take f passive , true( M ∗ , meas10 -µ ) = f passive , obs( M ∗ , meas), the exact solution for κ = 0. This results in a change to the definition of κ so that it no longer corresponds to an exact offset for active galaxy stellar masses; however, the correlation between κ and the average offset implied for active galaxies should still be high. With this approach, it is straightforward to calculate the measured stellar mass function:  \nSolving this equation is possible, but using the result has disadvantages. Specifically, calculating f passive , true couples the values of ill-constrained regions of the stellar mass function to well-constrained regions through this recurrence relation. The values for α -the faint-end slope of the SMHM relation-and κ will then be constrained by the MCMC algorithm attempting to improve fits to the well-constrained regions of the stellar mass function, and will lose their physical significance.  \nφ meas( M ∗ , meas) = f passive , true( M ∗ , meas10 -µ ) φ true( M ∗ , meas10 -µ ) + φ true( M ∗ , meas10 -µ -κ )(1 -f passive , true( M ∗ , meas10 -µ -κ )) . (C5)  \nOr, in terms of the observed passive fraction,  \nφ meas( M ∗ , meas) = f passive , obs( M ∗ , meas) φ true( M ∗ , meas10 -µ ) + φ true( M ∗ , meas10 -µ -κ )(1 -f passive , obs( M ∗ , meas10 -κ )) . (C6)  \nThen, to correct for observational completeness, we multiply by the appropriate factor:  \nφ obs( M ∗ , meas) = c ( z ) φ meas( M ∗ , meas) .  \n(C7)  \nFinally, as surveys usually return the stellar mass function over a range of redshifts (e.g., z 1 < z < z 2), we must perform one last integration to determine the survey mass function:  \nφ survey( M ∗ , meas) = ∫ z 2 z 1 φ obs( M ∗ , meas) dVc ( z ) Vc ( z 2) -Vc ( z 1) , (C8)  \nwhere Vc ( z ) is the comoving volume out to redshift z .', 'C.2. Calculating the Observed Specific Star Formation Rates and Cosmic Star Formation Rate Density': "Eq. B1 gives the average star formation rate for halos as a function of halo mass and redshift. Eq. 3 gives the median stellar mass for halos as a function of halo mass and redshift. Ignoring scatter, the fiducial SSFR is simply:  \nSSFR fid( M ∗ ) = SFR ( Mh ( M ∗ )) M ∗ . (C9)  \nScatter in the observed stellar mass at fixed halo mass results in contributions from lower- and higher-mass halos, weighted by number density. If there is no correlation between stellar mass and star formation rate beyond the correlation between halo mass and star formation rate, then the specific star formation rate would be as follows:  \nSSFR ( M ∗ )nc = M -1 ∗ φ -1 true ( M ∗ ) ∫ ∞ -M ∗ φ fid( M ∗ + ∆ M ∗ ) SFR ( Mh ( M ∗ + ∆ M ∗ )) P ( ∆ M ∗ ) d ∆ M ∗ . (C10)  \nHowever, with a nonzero covariance ( ρ ) between SFR and stellar mass at fixed halo mass this equation is amended to become:  \nSSFR ( M ∗ )true = M -1 ∗ φ -1 true ( M ∗ ) ∫ ∞ -M ∗ φ fid( M ∗ + ∆ M ∗ )( SFR ( Mh ( M ∗ + ∆ M ∗ ))10 ρ ∆ ' M ∗ ) P ( ∆ M ∗ ) d ∆ M ∗ ∆ ' M ∗ =log 10 ( M ∗ + ∆ M ∗ M ∗ exp(0 . 5( ξ ln10) 2 ) ) . (C11)  \nThe formula for ∆ ' M ∗ is necessary because Mh ( M ∗ + ∆ M ∗ ) gives the halo mass at which the median stellar mass is M ∗ + ∆ M ∗ ; however, the covariance between SFR and SM is calculated relative to the average stellar mass at fixed halo mass. The exponential in the denominator (exp(0 . 5( ξ ln10) 2 )) is exactly the ratio of the average value to the median value of a log-normal distribution with standard deviation ξ .  \nThere are restrictions on the plausible range of covariance coefficient ρ . At low redshifts, a multiple regression analysis on group catalogue data from the SDSS (Wetzel et al. 2012) suggests that, at fixed halo mass, the average star formation rate rises by 0.23 dex per 1 dex increase in stellar mass. At high redshifts (e.g., z > 6), the covariance will be almost perfect (1 dex rise in SFR per 1 dex rise in stellar mass) because of the short times that galaxies have had to form stars. Between these two times, it is reasonable to assume that the covariance increases with increasing redshift, so we adopt the following formula to interpolate between z ∼ 0 and z /greaterorsimilar 6:  \nρ ( a ) = 1 + (4 ρ 0 . 5 -3 . 23) a + (2 . 46 -4 ρ 0 . 5) a 2 . (C12)  \nThis formula for ρ ( a ) gives ρ (0) = 1, ρ (1) = 0 . 23, and ρ (0 . 5) = ρ 0 . 5, where ρ 0 . 5 is a free parameter in our analysis.  \nBy comparison, the total (cosmic) star formation rate density is straightforward to calculate; it is simply the sum of all star formation in all halos:  \nCSFR true = ∫ ∞ -∞ SFR ( Mh ) dN d log Mh d log Mh . (C13)  \nTo calculate the measured quantities (as opposed to true quantities), we adopt the same conversion between measured stellar mass and true stellar mass as for the stellar mass function. However, we do not impose additional corrections for active galaxies; by definition, because recent star formation is more evident in such galaxies, it is somewhat easier to recover the star formation rate, as opposed to the overall stellar mass. Therefore, we only apply a uniform systematic correction and do so in the same way as for stellar masses:  \nSFR meas( Mh ) = SFR true10 µ . (C14)  \nTo calculate SSFR meas and CSFR meas, we would therefore replace SFR true with SFR meas in Eqs. C11 and C13. There is one more complicating factor; namely the burstiness / dustiness of high-redshift galaxies. To account for this, we add additional factors to calculate the observed SSFR and CSFR:  \nSSFR obs( M ∗ ) = 1 1 -b (1 -c ( z )) SSFR meas( M ∗ ) CSFR obs = [1 -(1 -b )(1 -c ( z ))] CSFR meas . (C15)", 'D. THE USE OF DOUBLE POWER LAWS TO FIT THE STELLAR MASS - HALO MASS RELATION': 'A double power law has been used to fit the stellar mass (or luminosity) to halo mass relation at least since Yang et al. (2003). While using it is simple and convenient, the quantity and quality of galaxy observations have improved dramatically since its first use. To test its continued usefulness, we have used abundance matching to extract the implied stellar mass - halo mass (SMHM) relation from our low-redshift stellar mass functions (Moustakas et al. 2013; Baldry et al. 2008). To account for the presence of scatter in observed stellar mass at fixed halo mass, we use a method based on Richardson-Lucy deconvolution to extract the median stellar mass as a function of halo mass. We assume a combined scatter of 0.2 dex.  \nThe deconvolved SMHM relation is shown in Fig. 20. The figure demonstrates that a double power-law can approximate the relation, but can be discrepant by up to 0.1 dex. Adding one more parameter improves the fit substantially, to a mean offset of less than 0.025 dex in our case. The discrepancies in the SMHM relation translate directly to discrepancies in the assumed stellar mass function. This may be seen in the right-hand panel of Fig. 20; the fits for the SMHM relations have been used to  \n<!-- image -->  \nFIG. 20.Left panel: The stellar mass-halo mass relation as obtained by deconvolution of observed data (Moustakas et al. 2013; Baldry et al. 2008) with an assumed scatter in observed stellar mass at fixed halo mass of 0.2 dex; the dashed and dashed-dotted lines show a double power law fit and the five-parameter fit used in this work. The double power law fit differs by up to 0.1 dex from the observationally-derived relation. The effects on the stellar mass function are shown in the right panel; namely, the fits shown in the left panel have been convolved with 0.2 dex scatter as well as with the halo mass function so as to compare to the original stellar mass functions. Again, differences of 0.1 dex are seen in the double power law case, which is well outside of the reported error bars.  \n<!-- image -->  \npopulate the halo mass function and reconvolved with 0.2 dex scatter to determine the resulting stellar mass functions. While the five-parameter fit that we use is largely within the error bars of the observational data, the double power law fit gives results that are again discrepant by up to 0.1 dex or 4 σ from the observed stellar mass function.  \nIt is tempting to ignore this issue: after all, the size of the systematic uncertainties in stellar mass is significantly larger ( ∼ 0.25 dex). However, if it is indeed the case that the deviations from the power law are driven by systematic uncertainties, then part of the value of the physical interpretation of the double power law is lost: the high-mass and low-mass slopes are now contaminated by using the incorrect parametrization. In our five-parameter fit, the high-mass dependence and low-mass slope are somewhat more insulated from unparametrized systematic errors in the transition region, leading to a more fair estimation of the high- and low-mass behavior of the SMHM relation.  \nNaturally, no fit will perfectly reproduce the original function. Special caution is necessary when using a MCMC method where a parametrization to tightly-constrained data is accompanied by additional systematic or nuisance parameters. This is because the MCMC method may try to recruit the additional parameters to improve the fit, rather than for their intended purpose (i.e., to cover the range of allowed uncertainties). As shown in Fig. 20, our chosen parametrization can reproduce the local SMF to within 0.02 dex in the vast majority of cases. Hence, to avoid overfitting bias, we reduce the magnitude of errors by 0.02 dex at each data point before calculating the χ 2 value used with our MCMC method. At z < 0 . 2, this has the impact of decoupling the systematic uncertainty parameters from the detailed shape of the SMF; at higher redshifts, there is no effect, as the typical size scale of the errors for observational data is already 0.1-0.2 dex.', 'E. THE REDSHIFT EVOLUTION OF THE STELLAR MASS - HALO MASS RELATION FROM ABUNDANCE MATCHING': "As a check of the accuracy of our analysis pipeline, it is instructive to compare with a direct abundance matching approach. In Fig. 21, we show results for both direct abundance matching (left panel) and abundance matching including an assumed intrinsic ( ξ ) and measurement ( σ ) scatter that match our assumed fiducial defaults. Notably, several features present in our main analysis are evident in this simpler approach:  \n- 1. The SMHM relation for low-mass halos ( Mh < 10 11 M /circledot ) at z = 0 deviates from a power law.\n- 2. The halo mass at which the integrated star formation is most efficient rises by 0.5 dex from 10 12 M /circledot at z = 0 to 10 12 . 5 M /circledot at z = 2.\n- 3. Above z = 2 it is difficult to locate the halo mass at which the peak integrated efficiency occurs because surveys at high redshifts do not have sufficient volume to probe the exponential tail of the stellar mass function (if one exists). If the UV luminosity function is any indication, however, the z = 8 results suggest that the peak may occur at lower masses above z = 2.\n- 4. There appears to be a turnaround in the SMHM relation near z ∼ 2. The integrated efficiency of Mh ∼ 10 11 . 5 M /circledot halos falls from z = 0 to z = 2, but then increases at higher redshifts.  \nAt high redshifts, it is difficult to tell what, if any, evolution there is in the shape of the SMHM relation. However, it is clear that either the efficiency or the characteristic halo mass scale is evolving. If we constrain the shape parameters in our fitting formula  \n<!-- image -->  \nFIG. 21.Left panel: the stellar mass to halo mass ratio inferred via direct abundance matching of stellar mass functions (Table 3) to halo mass functions from z = 0 to z = 8, with zero scatter. Right panel: same, with deconvolution of 0.16 dex intrinsic scatter (in stellar mass at fixed halo mass) and (0 . 07 + 0 . 04 z ) dex of measurement scatter (in measured stellar mass at fixed true stellar mass) applied. Note that deconvolution applied to noisy data (e.g., SMFs) amplifies the noise, which is why features in the left panel appear smoother than those in the right panel.  \n<!-- image -->  \nTABLE 6 EMPIRICAL FIT TO COSMIC STAR FORMATION RATES  \nNOTE. - Fits according to Eq. F1 for the collected cosmic star formation rate data in Table 4 ('This Work') and the data in Hopkins & Beacom (2006), all corrected to a Chabrier (2003) IMF. Note that these are fits to observational data, not fits to our best-fit model. The quality of the fits may be seen in Fig. 2. The constant C has units of M /circledot yr -1 comoving Mpc -3 ; all other constants are unitless. Systematic errors are estimated in Table 7.  \nTABLE 7 ERROR ESTIMATES FOR THE COSMIC STAR FORMATION RATE  \nNOTE. - 'Inter-publication' refers to the standard deviation in dex between published data points (from references in Table 4) in redshift bins after correction for the evolution of the cosmic star formation rate (Eq. F1); the error in the error estimate for each redshift bin is of order 25%. 'This Work - HB06' refers to the difference in dex between the best fit (Eq. F1) to data points in Table 4 and the best fit to data points in Hopkins & Beacom (2006) over the same redshift bins; this is a qualitative estimate of how much assumptions about calculating the cosmic star formation rate have changed over the past decade.  \nfor the SMHM relation ( α,δ, and γ ) to the z = 0 values, we find that high redshifts strongly prefer a halo mass scale that scales as a power law with redshift, but have comparatively little implied evolution in the efficiency at high redshifts. For that reason, when we extrapolate to high redshifts unconstrained by data in our main analysis ( z > 8 . 5), we assume that /epsilon1, α,δ, and γ do not evolve beyond z = 8 . 5, but that the characteristic halo mass M 1 continues to evolve according to the redshift scaling we have adopted.", 'F. CALCULATINGUNCERTAINTIESIN COSMIC AND SPECIFIC STAR FORMATION RATES': "For cosmic star formation rates (SFRs), we selected papers from the astro-ph arXiv that  \n- 1. Were posted within the past 6 years (2006-2012).\n- 2. Contained the word 'cosmic' and either the word 'SFR' or the words 'star' and 'formation' in the title or abstract.\n- 3. Contained an estimate of the cosmic star formation rate from observations.\n- 4. Were not superseded by later publications.\n- 5. Contained a reference to the initial mass function (IMF) used.  \nWe have in addition excluded constraints from long gamma-ray burst SFR estimates for this work (e.g., Kistler et al. 2009); we feel that not enough data points are available at this time to independently calibrate the evolution of known systematic biases (e.g., higher GRB rates in low-luminosity galaxies) with redshift. The selected papers are summarized in Table 4. From these,  \n<!-- image -->  \nFIG. 22.- Plots of the specific star formation rate data assembled in this study (see Table 5) for galaxies of stellar mass 10 9 . 5 M /circledot (top) and 10 10 . 5 M /circledot (bottom). Double power-law fits (Eq. F1) to the data sets are taken from Table 8. Error bars are not shown since authors often dramatically underestimate the true magnitude of their systematic errors (see §3.3).  \n<!-- image -->  \nTABLE 8 EMPIRICAL FITS FOR SPECIFIC STAR FORMATION RATES  \nNOTE. - Fits according to Eq. F1 for the collected specific star formation rate data in Table 5, all corrected to a Chabrier (2003) IMF. The quality of the fits may be seen in Fig. 22. Note that these are fits to observational data, not fits to our best-fit model. Asterisks denote that the parameter in question ( z 0) was consistent with 1.0 when allowed to vary; hence, all fits were constrained to have z 0 = 1 so as to allow better comparison between the evolution of the remaining fit parameters with stellar mass. The constant C has units of yr -1 ; all other constants are unitless. Systematic errors are estimated in Table 9.  \nTABLE 9 ERROR ESTIMATES FOR SPECIFIC STAR FORMATION RATES  \nNOTE. - The standard deviation in dex between published data points (from references in Table 5) in redshift bins after correction for the evolution of the specific star formation rate (Eq. F1). For reasons mentioned in §3.3, the average does not include the estimates in the z = 3 -8 redshift range. The error in the error estimate for each bin is of order 45%; the error on the average error is of order 12%. Errors are not estimated when there are fewer than four references in a given stellar mass and redshift range.  \nwe calculate the inter-publication variance in five redshift bins with 16 data points each. To take into account variation of the cosmic SFR in bins, we subtract a double-power law fit to the data points in Fig. 2:  \nCSFR ( z ) = C 10 A ( z -z 0) + 10 B ( z -z 0) , (F1)  \nwhere A , B , C , and z 0 are constants given in Table 6; adding more parameters did not make the quality of the fit substantially better. Our results are given in Table 7; we find the average systematic error to range from 0.13 dex at z = 0 to 0.27 dex for z > 3. For SSFRs we conducted a very similar literature search. In this case, the selection criteria on astro-ph were papers that were posted within the past 6 years (2006-2012); contained the word 'specific' and either the word 'SFR' or the words 'star' and 'formation' in the title or abstract; contained an estimate of the specific star formation rate as a function of stellar mass, both as calculated from observations; were not superseded by later publications; contained a reference to the initial mass function (IMF) used; and did not limit themselves to a specific class of galaxies (sub-mm, ULIRGs, satellites, star-forming only, etc.). The selected papers are summarized in Table 5.  \nTo estimate uncertainties, we interpolate data for each publication's tabulated redshift ranges to obtain (if possible) estimates of the SSFR at four stellar masses (10 9 , 10 9 . 5 , 10 10 and 10 10 . 5 M /circledot ). We divide the SSFR data into the same five redshift bins as for the cosmic SFR data and estimate the systematic errors from the inter-publication variance. To exclude the effects of evolution in the SSFR over a redshift bin, we again subtract a double-power-law fit (Eq. F1) with parameters determined separately for each  \na  \nFIG. 23.Left: A comparison between the Tinker et al. (2008) mass function and two simulations with different initial conditions ( ZA : Zel'dovich Approximation; 2LPT : second-order Lagrangian Perturbation Theory). Specifically, the plot shows the ratio of the cumulative central halo mass function to 10 11 . 5 M /circledot ( Φ (10 11 . 5 M /circledot )) in the simulations to that for the Tinker MF. Small mass-dependent corrections (Eq. G1) were made to renormalize the Consuelo mass function to match the Tinker MF at z = 0, but no such corrections were necessary for Bolshoi. The Tinker MF, which was calibrated using ZA simulations, agrees remarkably well with Bolshoi over a wide range of redshifts; however, it is noticeably different from the Consuelo simulation for z > 3. As several authors (Crocce et al. 2006; Jenkins 2010) have shown improved convergence for 2LPT over ZA initial conditions, we use the 2LPT normalization as the basis for our corrective fit (Eq. G2). Right: the mass-dependence of the correction to the Tinker cumulative mass function from the Consuelo simulation. Solid lines show direct results from the simulation; dashed lines show the fit in Eq. G3. The dot-dashed line shows the fit in Eq. G3 extrapolated down to Φ (10 10 M /circledot ).  \n<!-- image -->  \n<!-- image -->  \na  \nstellar mass in Table 8; the quality of the fit is shown in Fig. 22. It is worth remarking that Weinmann et al. (2011) conduct a similar analysis of SSFRs at M ∗ = 10 9 . 5 M /circledot and find a 'plateau' at high redshifts; while the data are consistent with a plateau, they are also consistent with not having a plateau: indeed, any range of shallow slopes at high redshifts would be consistent with the data collected here.  \nOur results for the inter-publication variance as a function of stellar mass and redshift are shown in Table 9. Because only UV surveys can probe low stellar masses at high redshifts ( z > 3), the published SSFRs for low stellar masses at high redshift are more in agreement (0.22 dex scatter) than for low stellar masses at low redshift (0.32 dex scatter)! This strongly suggests that the inter-publication variance in these bins does not represent the true systematic uncertainties. Excluding high-redshift bins ( z > 3), there are no obvious systematic trends in the inter-publication variance, either with redshift or stellar mass. Indeed, a one-way analysis of variance (ANOVA) test delivers no strong evidence that the inter-publication variances differ across bins ( p = 0 . 6). This conclusion is not affected by merging bins by stellar mass or by redshift, and the size of the effect is too small to be measurable without adding more data. We adopt the fiducial model of constant errors across stellar masses and redshifts, where the error is determined by merging all redshift and stellar mass bins except for the high-redshift bins ( z > 3) which are known underestimates, we would estimate the average systematic error in SSFRs to be 0.28 dex.  \nFinally, we mention a commonly-overlooked systematic error in optical SFRs, which comes from the bias introduced by separately fitting the dust model for each galaxy. As the dust model for each galaxy can never be derived exactly, individual galaxy SFRs will have some scatter relative to the true galaxy SFRs. One way to estimate the shape of this scatter is to compare two different SFR estimates with different assumptions; for example, Salim et al. (2007) finds log-normal scatter consistent with intrinsic scatter of 0.3 dex. In calibrating one estimator to another (or to observations), authors typically make sure that the median offsets are zero; however, they often fail to notice that the expected mean offset is nonzero because the log-normal distribution is skewed in linear space. Specifically, a log-normal distribution with median value 1 and scatter σ (in dex) has an expected mean value of  \n〈 x 〉 = exp(0 . 5( σ ln10) 2 ) . (F2)  \nFor a log-normal distribution with 0.3 dex scatter, this corresponds to an overall bias of +0.1 dex. This correction does not apply to UV studies that apply a median dust correction to all galaxies or to those that already correct for this bias.", 'G. CORRECTIONSTO AND UNCERTAINTIESIN THE TINKER MASS FUNCTION': "The halo mass function for distinct halos presented by Tinker et al. (2008) was calibrated for a wide array of simulations from z = 0 to z = 2 . 5. The present work requires a mass function from z = 0 to z > 8 including satellite halos, so some additional corrections must be made. In addition, as new studies (Reddick et al. 2012) suggest that abundance matching to halo properties at their peak value before stripping offers a better match to z = 0 clustering data, calibration for the quantities of interest using M peak must also be applied.  \nTo test the high-redshift calibration of the mass function in Tinker et al. (2008), we integrate the mass function in two simulations (Bolshoi and Consuelo) to 10 11 . 5 M /circledot and compare to the Tinker mass function for the appropriate cosmology. A massdependent correction is necessary to normalize the z = 0 mass function for Consuelo to the Tinker mass function to account for slight incompleteness in the Consuelo mass function at this level:  \nlog 10 ( φ Consuelo( M ) φ Tinker( M ) ) = -0 . 0008 -0 . 042 ( M 10 11 . 51 M /circledot ) -0 . 69 (G1)  \n<!-- image -->  \nFIG. 24.- Best fits to the evolution of the satellite fractions (Eqs. G7-G10) in Bolshoi for satellite accretion masses (left panel) and satellite peak masses (right panel), according to the simplified model of Eq. G6. The colored lines correspond to the fitted satellite fractions; the black dashed lines correspond to the satellite fractions in Bolshoi.  \n<!-- image -->  \nThe net effect of this correction on Φ (10 11 . 5 M /circledot ) is an increase of 0.03 dex (7%) compared to the raw simulation result for Consuelo. No correction is necessary for the Bolshoi simulation. As shown in Fig. 23, the redshift fit given in Tinker et al. (2008) matches Bolshoi to within 0.02 dex up to a redshift of 5-6. However, Consuelo diverges from the Tinker mass function by 0.06 dex (15%) already at z = 3. This discrepancy is due to the different initial conditions calculations used in the two simulations: Consuelo uses second-order Lagrangian Perturbation Theory (2LPT), whereas Bolshoi and the Tinker MF use the Zel'dovich approximation (1LPT). The Zel'dovich approximation has been previously found to underestimate the early nonlinear collapse of overdense regions (Crocce et al. 2006; Jenkins 2010), which is reconfirmed by the comparison in Fig. 23. For that reason, we prefer to use the correction from Consuelo, parametrized as follows:  \nlog 10 ( Φ true(10 11 . 5 M /circledot ) Φ Tinker(10 11 . 5 M /circledot ) ) = 0 . 144 1 + exp[14 . 79( a -0 . 213)] . (G2)  \nAs shown in Fig. 23, this fit matches the correction from the Consuelo simulation to within 1%. However, as also shown in Fig. 23, there is a slight mass dependence to the correction; the magnitude of the correction is larger for more massive halos, as their nonlinear collapse begins earlier. We parametrize the mass dependence as follows:  \nlog 10 ( Φ true( M ) Φ Tinker( M ) ) = log 10 ( Φ true(10 11 . 5 M /circledot ) Φ Tinker(10 11 . 5 M /circledot ) )( M 10 11 . 5 M /circledot ) 0 . 5 1 + exp(6 . 5 a ) . (G3)  \nThis function is well-behaved, even when extrapolated down to low masses, as shown in Fig. 23.  \nParametrizing the corrections as we have done introduces several uncertainties. The use of only a single simulation to calibrate the mass function introduces sample variance uncertainties, but these are in fact small for a box as large as Consuelo (420 Mpc h -1 ); at low redshifts, the median offset from the Tinker mass function is 0.2% (see Eq. G1), and even at redshift 6, the Poisson and sample variance errors combined are on the order of 1% (Trenti & Stiavelli 2008) for Φ (10 11 . 5 M /circledot ). Much larger uncertainties come from the use of a mass-dependent incompleteness correction to the Consuelo mass function, because uncertainties in this correction are degenerate with the mass dependence of the normalization correction in Eq. G3. If the incompleteness correction needed to be larger at high redshift, the evidence for the mass-dependence in the normalization correction would weaken. We can take this as an estimate of the additional systematic uncertainties on the high-redshift mass function. In terms of the impact on the halos in which high-redshift galaxies ( z =7 -8) are expected to reside, the number density of such objects would be underestimated by at most 0.05 dex (12%), according to Fig. 23, if there were no mass-dependence in the normalization correction. This is tiny in comparison to the uncertainties in the stellar masses and star formation rates at high redshift. Indeed, even if the real universe corresponded exactly to the Tinker mass function, the total error introduced by using the corrections in this section (0.1 dex, or 26%) would again be dwarfed by uncertainties in measuring galaxy stellar content.  \nIn terms of the satellite fraction, we determine fits using a simple, physically-motivated model. It is an oft-observed feature of dark matter substructure that the number of satellite halos at a given mass scales approximately with the mass of the host the number density of satellites in hosts of mass M c is proportional to M c φ c( M c), where φ c ( M ) is the mass function for central halos. Another feature of dark matter substructure is that the number of satellites at a given satellite-to-host mass ratio is largely self-similar across host halo masses; this requires that the proportionality must scale inversely with the satellite mass. Therefore, we can write down a simple formula for the satellite mass function, φ s :  \nφ s ( M ) ∼ C ( a ) M ∫ ∞ M M c φ c( M c) d log 10 M c (G4)  \nwhere the proportionality factor, C ( a ), is a function of scale factor only. The low-mass limit of this function takes on an especially  \nwhich can be simplified to  \nφ s ( M ) φ c ( M ) ∼ C ( a ) log 10 ( M cutoff( a ) M ) (G6)  \nThus, the satellite fraction grows only logarithmically in the low-mass limit.  \nHowever, the picture in Eq. G4 is not entirely complete; because they have later assembly times, massive clusters are more likely to be undergoing major mergers than smaller halos. However, these major mergers also become absorbed more quickly, as the effects of dynamical friction are more pronounced in massive halos. These two (opposite) effects combine to give a satellite fraction slightly larger at the high-mass end than Eq. G4 would predict; moreover, this excess evolves nontrivially with redshift. For that reason, we use the slightly more flexible form of Eq. G6 and allow great freedom in the parametrization for both C ( a ) and M cutoff( a ).  \nFitting the satellite fraction to a simulation is made especially tricky by the fact that simulations are incomplete for satellites at higher masses than for central halos, even before selecting satellites on M peak or M acc. For that reason, Consuelo is unsuitable for determining the fraction of all but the most massive satellites ( M ∼ 10 13 M /circledot ), which provides an especially poor handle on constraining the redshift evolution of the fit. Bolshoi fares better, being complete for satellites down to 10 11 M /circledot ; yet, given its use of ZA initial conditions, it may not seem especially suitable for high-redshift calibration, either!  \nYet, there are two arguments in its favor. First, the satellite fraction is low at high redshift, so even large relative errors in the satellite fraction will translate into small errors on the overall mass fraction. Second, the influence of the initial conditions on the satellite fraction are largely through the effects of assembly bias; however, recent work has shown that satellite richness (when tallied above a fixed cut in accretion or peak mass) is almost uncorrelated with formation time (Wu et al., in prep.). The impact of the initial conditions on the satellite fraction is thus minor compared to the effect of the mass function normalization.  \nTherefore, using the Bolshoi simulation, we fit the satellite fraction from z = 0 to z = 10; when satellite masses are determined by M acc, we find:  \nlog 10( C ( a )) = -2 . 69 + 11 . 68 a -28 . 88 a 2 + 29 . 33 a 3 -10 . 56 a 4 (G7)  \nlog 10 ( M cutoff( a )) = 11 . 34 + 8 . 34 a -0 . 36 a 2 -5 . 08 a 3 + 0 . 75 a 4 (G8)  \nand find the following parametrization for the satellite fraction when satellite masses are determined by M peak:  \nlog 10 ( C ( a )) = -1 . 91 + 6 . 23 a -15 . 07 a 2 + 15 . 02 a 3 -5 . 29 a 4 (G9)  \nlog 10 ( M cutoff( a )) = 10 . 66 + 15 . 93 a -21 . 39 a 2 + 18 . 20 a 3 -8 . 21 a 4 (G10)  \nBoth of these fits in comparison to satellite fractions in Bolshoi are shown in Fig. 24; they are accurate to within the estimated uncertainties in the mass function (5%, Tinker et al. 2008).", 'H. HALO MASS ACCRETION RATES': 'The most commonly used fitting formula for halo mass accretion rates is that proposed by Wechsler et al. (2002), which suggests that halo mass growth is an exponential function in redshift ( M ( z ) ∝ exp( -cz )). As simulations improved, deviations from this simple formula were found, which led authors to propose more complicated formulas; e.g., M ( z ) ∝ (1 + z ) b exp( -cz ) in McBride et al. (2009). However, little attention has been given to the high-redshift behavior of these formulas, with the predictable result that they fail to accurately capture median mass accretion histories for halos as recently as z ∼ 3 (see Fig. 25, left panel). 6  \nThis limitation a problem for this study. To motivate our choice of improved fitting formulas, we note an issue with all proposed fitting formulas to date, which is their exponential dependence at high redshifts. That is to say, the mass ratio of two halos M 1 and M 2 will asymptote to  \nM 1( z ) M 2( z ) → M 1(0) M 2(0) exp(( c 2 -c 1) z ) (H1)  \nThis is problematic because more massive halos usually have steeper exponential dependencies (i.e., M 1 < M 2 implies c 1 < c 2). If M 1( z = 0) < M 2( z = 0), these formulas then predict that eventually the progenitor histories will cross and M 1( z ) > M 2( z ). Using the McBride et al. (2009) formula, for example, with the MultiDark simulation, we find that it predicts the progenitors of 10 15 M /circledot halos to be less massive than the progenitors of 10 14 . 75 M /circledot halos as early as z = 5. However, this is not the case: in simulations, the hierarchy of average or median progenitor masses is maintained, even though the distance between them may decrease (see, for example, the right panel of Fig. 25).  \nWe therefore parametrize the mass accretion history of a single halo mass, and then parametrize the ratios of progenitor masses for other halo masses, which entirely avoids the problem of crossing progenitor mass histories. The mass bin we choose  \ncision is not as necessary for cases when individual stochasticity in merger events results in higher deviations than those shown in the bottom-left panel of Fig. 25.  \nappealing form; as φ c ( M c) ∼ φ 0 M -1 c for Mc is less than the exponential cutoff scale, M cutoff, we find:  \nφ s ( M ) φ c ( M ) ∼ C ( a ) M ∫ M cutoff M φ 0 d log 10 Mc φ 0 M -1 (G5)  \nFIG. 25.- Fits to median mass accretion histories (Mvir) for central halos in Bolshoi, MultiDark, and Consuelo simulations. Common functional forms for the mass accretion rates, such as the exponential fit of Wechsler et al. (2002) or the exponential plus power-law fit of McBride et al. (2009) fail to accurately capture mass accretion histories beyond z = 2. Moreover, at high redshifts, these common functional forms imply unphysical behavior-i.e., that progenitors of 10 15 M /circledot halos are on average less massive than progenitors of 10 13 M /circledot halos (see text), which is not at all the case. The left panel shows mass accretion histories for 10 13 . 15 to 10 13 . 40 M /circledot halos in Bolshoi and Consuelo, along with residuals (linear scale). Bolshoi and Consuelo have different cosmologies, so it is expected that the mass accretion rates of halos differ slightly. The right panel shows progenitor mass ratios (compared to progenitors of 10 13 M /circledot halos) for halos of mass 10 11 to 10 15 M /circledot at z = 0 for the Bolshoi, MultiDark, and Consuelo simulations, along with residuals. In all cases, our fits to the mass accretion histories of halos for the Bolshoi cosmology are accurate to 5% on average; the only exception is for the highest mass bins in the Bolshoi and Consuelo simulations, which contains only a few halos.  \n<!-- image -->  \n<!-- image -->  \na  \n(10 13 . 15 -10 13 . 40 M /circledot ) is well-resolved in all three dark matter simulations at our disposal, and sufficient statistics are available even in the smallest simulation (Bolshoi) to be able to reduce scatter to an acceptable level. Our resulting fit is:  \nM ( M 0 , z ) = M 13( z )10 f ( M 0 , z ) (H2)  \nM 13( z ) = 10 13 . 276 (1 + z ) 3 . 00 (1 + z 2 ) -6 . 11 exp( -0 . 503 z ) M /circledot (H3)  \nf ( M 0 , z ) = log 10 ( M 0 M 13(0) ) g ( M 0 , 1) g ( M 0 , 1 1 + z ) (H4)  \ng ( M 0 , a ) = 1 + exp( -4 . 651( a -a 0( M 0)) (H5)  \na 0( M 0) = 0 . 205 -log 10 [ ( 10 9 . 649 M /circledot M 0 ) 0 . 18 + 1 ] (H6)  \nwhere M ( M 0 , z ) gives the median virial mass (Bryan & Norman 1998) for progenitors of halos with mass M 0 at z = 0. Despite its complexity, this seven-parameter fit (three parameters for M 13( z ), one parameter for g ( M 0 , a ), and three parameters for a 0( M 0)) provides an excellent fit to mass accretion histories for a wide range of halo masses and redshifts, as shown in Fig. 25.', 'I. SUBHALOMERGER/DISRUPTIONRATES': "Halo mergers are important for our study because they impact the rate at which stellar mass is transfered from low-mass galaxies to high-mass galaxies as well as the rate at which stars are ejected into the intracluster light (ICL). While there have been many studies of the rate at which halos become subhalos (see, e.g., Fakhouri et al. 2010) in dark matter simulations, these studies are less relevant to our work because it is expected that satellite galaxies should last at least as long as the corresponding subhalo that hosts them (e.g. Yang et al. 2012; Reddick et al. 2012). Studying the subhalo destruction rate is more difficult because of issues with halo finding completeness when subhalos pass close to the centers of their hosts (Behroozi et al. 2013b,c), and relatively few authors have attempted to correct for this issue in their work (c.f., Wetzel & White 2010).  \nHowever, with the advantage of a phase-space halo finder and a merger tree algorithm that corrects for such inconsistencies, we are in a unique position to provide a calibration of the merger rate (or disruption rate) of subhalos into host halos all the way back to the earliest redshifts in simulations. We find that the functional form proposed in Fakhouri et al. (2010) is an excellent fit to the subhalo disruption rate; however, the parameter values and redshift evolution are different due to the different merger definition, cosmology, and mass definition used.  \nLetting Mh be the mass of the host halo and θ be the ratio of subhalo peak mass to host mass, we find that the disruption rate (expressed in terms of the decrement in number of subhalos per unit host halo per unit redshift per log interval in subhalo mass  \nFig. 26 shows that scaling in both mass and redshift of the subhalo disruption rate has been successfully captured in Eqs. I1-I3. The remaining deviations are due to Poisson noise in the disruption rates for halo masses where Bolshoi has few halos (especially in the highest mass bin) and stochasticity in halo mass estimations for low-mass halos. Specifically, as the peak mass captures the highest mass ever recorded in a halo's mass accretion history, greater variance in the halo mass estimator translates to artificially increased peak masses for low-mass halos. This manifests as a slight upturn in the subhalo disruption rate close to the halo mass resolution limit of the simulation.  \n<!-- image -->  \nFIG. 26.- The rate of subhalo disruption in host halos in the Bolshoi simulation with mass and redshift dependencies ( X ( M , z )) scaled out. Results are presented in terms of the number of subhalos disrupted per host halo per unit redshift per log interval in satellite mass ratio. The subhalo merger mass ratio ( θ ) in these plots is the ratio of the subhalo peak mass (i.e., the highest mass any progenitor of the subhalo ever had) to the host halo mass. The left panel shows the scaled subhalo disruption rate at z ∼ 0 over a range of host halo masses, demonstrating excellent convergence. The red dashed line in the left panel represents the average of scaled subhalo disruption rates as a function of subhalo mass ratio θ , excluding halos for which subhalos of the given mass ratio would be below the mass resolution limit of the simulation (10 10 M /circledot ). The right panel shows the mass-averaged scaled subhalo disruption rate as a function of redshift, as well as the fit we adopt in this work. The fit shows excellent conformance to the simulation data, except for very low subhalo mass ratios, where stochasticity in evaluating halo masses at low particle numbers artificially inflates peak halo masses (see text).  \n<!-- image -->  \nratio) is:  \n-d 2 N dzd log 10 θ ( M , θ, z ) = A ( M , (1 + z ) -1 ) θ -1 exp( b ( M , θ )) (I1)  \nA ( M , a ) = 0 . 0316 ( 1 + 0 . 20ln( a ) 1 + a )( M 10 12 M /circledot ) 0 . 03 + 0 . 05 a (I2)  \nb ( M , θ ) = [ 10 1 . 929 ( M 10 12 M /circledot ) 0 . 11 θ ] 0 . 2586 (I3)  \nIf θ is instead the ratio of subhalo accretion mass to host halo mass, we find instead a slightly simpler dependence:  \nA acc( M , a ) = 0 . 0139 ( 1 + 0 . 25ln( a ) 1 + a )( M 10 12 M /circledot ) 0 . 052 + 0 . 072 a (I4)  \n[ ] We show the quality of the fits for the peak subhalo mass in Fig. 26. For ease of comparing simulation data to the fit, we have scaled out host mass and redshift dependencies according to Eqs. I1-I3, so that all plotted lines are directly comparable to the subhalo disruption rate for M = 10 12 M /circledot at z = 0.  \nb acc( M , θ ) = 10 3 . 091 θ 0 . 1965 (I5)", 'J. FITS TO STAR FORMATION RATES AND HISTORIES UNDER DIFFERENT MODELING CONSTRAINTS': "We consider four models with alternate sets of assumptions to the main analysis in this paper:  \n- 1. A model ('No Systematics') that excludes all nuisance parameters related to errors in converting luminosities into stellar masses (i.e., µ = κ = 0, η cal = ∞ ) and to incompleteness ( A = 0).\n- 2. A model ('No z > 6 Constraints') that excludes all data constraints at high redshift ( z > 6), including data in Bouwens et al. (2011); Bradley et al. (2012); Schaerer & de Barros (2010); McLure et al. (2011), and the data points above z = 6 in Bouwens et al. (2012).\n- 3. A model ('No CSFR Constraints') that excludes all constraints from observed cosmic star formation rates.\n- 4. A model ('No SFR Constraints') that excludes all constraints from observed star formation rates (both CSFRs and SSFRs) as well as excluding all systematic nuisance parameters as in the first alternate model ('No Systematics').  \n<!-- image -->  \n<!-- image -->  \n<!-- image -->  \nFIG. 27.- Comparisons to observational data for the best fits of alternate models. Top left : comparisons to observed stellar mass functions for the model without systematic uncertainties (i.e., µ = 0, κ = 0, η cal = ∞ , no incompleteness). Top right : comparisons to stellar mass functions for a model excluding all constraints at z > 6 (i.e., fromBouwens et al. (2011); Bradley et al. (2012); Bouwens et al. (2012); Schaerer & de Barros (2010); McLure et al. (2011)). Bottom left : comparisons to observations of the cosmic star formation rate for the two models already mentioned, along with a model that excludes all cosmic star formation rate constraints ('No CSFR'; using SSFRs only), and a model that excludes all star formation rate constraints and systematic effects ('No SFR'). Bottom right : comparisons with observed SSFRs for all five models.  \n<!-- image -->  \nComparisons with observed data are presented in Fig. 27. For ease of comparison across mass ranges, we show comparisons between derived baryon conversion efficiencies and stellar mass histories for the main model in this paper and the alternate models in Figs. 28 and 29.  \nFirst, it may be verified that the model including no systematics succeeds remarkably well in matching all three kinds of observations; while it is potentially somewhat high at low redshifts for the CSFR, this is not unexpected on account of issues with the galaxy / ICL definition (§2.4). The No Systematics model has much tighter error bars on the stellar mass to halo mass ratios (not shown) than the full model at all redshifts, which suggests that systematics remain the single largest source of errors in determining stellar masses. Fig. 28 is comparable to Fig. 5 for the No Systematics model, except that at high redshifts and low halo masses, there is an upturn in the SFR that is not present in the full model.  \nThe model excluding z > 6 data gives almost identical results as the full model for z ≤ 6, which is encouraging. The 'predictions' for high-redshift stellar mass functions are somewhat lower than in the best-fitting model, which results in a higher SFR (i.e., rate of stellar mass growth) over the period from z ∼ 7 to z ∼ 6.  \nThe model excluding all cosmic star formation constraints is remarkably similar to the full model in our analysis. This is as expected: the specific star formation rate in combination with the stellar mass function gives a constraint on the total amount of star formation. Yet, the CSFR data has tighter error bars than the SSFR in combination with the SMF, at least for z < 1 (see §3.2 and 3.3), so it provides slightly better constraints at those redshifts.  \nExcluding all star formation rate data results in a good fit to the CSFR at high redshifts, where growth in stellar mass is rapidi.e., where stellar mass functions at successive redshifts are significantly different even despite large errors. At lower redshifts, the predicted cosmic star formation rate for the best-fitting model is discrepant from observations. At these redshifts, the growth in stellar mass is poorly constrained because the intrinsic stellar mass functions are changing comparatively little, especially from  \nz  \n<!-- image -->  \n<!-- image -->  \nz  \n<!-- image -->  \nFIG. 28.- Derived baryon conversion efficiencies and stellar mass histories for three alternate models (excluding all systematic uncertainties, excluding constraints at z > 6, and excluding constraints from cosmic star formation rates); the shaded bands show the one-sigma distributions from the main model in this paper, and the solid lines show results from the alternate models.  \n<!-- image -->  \n<!-- image -->  \n<!-- image -->  \n<!-- image -->  \nFIG. 29.- Derived baryon conversion efficiencies and stellar mass histories for one alternate model (excluding all systematic uncertainties and constraints from star formation rates); the shaded bands show the one-sigma distributions from the main model in this paper, and the solid lines show results from the alternate model.  \n<!-- image -->  \nz = 0 . 5 to z = 0. A small difference in the choice of stellar mass function evolution at low redshifts therefore has a much larger impact on the relative error in star formation rates. Thus, although the stellar mass to halo mass ratios are similar to the No Systematics model (and the abundance matching model in Appendix E), the star formation rates derived using this model are unreliable."}