EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 615 (July 2018) A&A, 615 (2018) A146 Full HTML
Free Access
Issue
A&A
Volume 615, July 2018
Article Number A146
Number of page(s) 15
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/201832821
Published online 30 July 2018

© ESO 2018

1 Introduction

It has been observed that there is a strong correlation between the stellar mass (M⋆) and the star formation rate (SFR) for the majority of star forming galaxies (SFG; e.g. Brinchmann et al. 2004; Noeske et al. 2007; Elbaz et al. 2007; Speagle et al. 2014; Tomczak et al. 2016). This has become known as the main sequence (MS) of star forming galaxies.

The MS is notable due to its tight correlation; the scatter of the SFR-M⋆ relation has been found to be approximately 0.3 dex. This low scatter has been observed to exist for approximately the last 10 Gyr and has been found to be independent of M⋆ (Whitaker et al. 2012; Speagle et al. 2014; Tomczak et al. 2016). The consistency of the scatter is understood to be a result of every galaxy having star formation (SF) regulated by similar quasi-static processes (Lee et al. 2015), while the scatter arises as a result of minor fluctuations of the flow of material into galaxies (Tacchella et al. 2016; Mitra et al. 2017). Galaxies move to the upper MS as the gas is compacted in the centre of the galaxy which triggers SF. As the central gas is depleted, but before the galaxy moves above the MS, the SFR reduces and the galaxy falls to the lower MS. New gas then falls into the galaxy, replenishing the central gas reservoir before the galaxy becomes quiescent. This cycle repeats until the galaxy’s gasreplenishment time is longer than the depletion time and the galaxy becomes quiescent (Tacchella et al. 2016).

Two different schools of thought exist over the shape the MS takes. Some studies find the relation between SFR and stellar mass is a simple power law (e.g. Speagle et al. 2014), i.e. the log of the SFR increases with the log of the stellar mass as (1)

where α is the slope and β is the normalisation. However, other studies find there is a high-mass turn-over at log(M⋆/M⊙) ≈ 10.5 with the slope of the MS being shallower above the turn-over than below (e.g. Lee et al. 2015; Tomczak et al. 2016). Recent studies have shown the two different forms of the MS may be a result of how the SFG population is separated from the quiescent galaxies (QG). Johnston et al. (2015) have shown that using a SFG-QG cut that is stricter in selecting the SFG population, in this case a 4000 Å break index value that is lower, will result in a straight MS, while a cut that is less strict forms a MS with a turn-over. The less strict cut leaves in high-mass, low SFR objects, which lower the mean SFR at high mass.

Regardless of the form of the MS, the normalisation is found to increase with redshift (e.g. Speagle et al. 2014; Johnston et al. 2015; Schreiber et al. 2015; Tomczak et al. 2016). This increasing normalisation is not surprising. As redshift increases, the fraction of cold gas available for star formation increases (Tacconi et al. 2010; Dunne et al. 2011; Genzel et al. 2015; Scoville et al. 2016). Thus, with more gas available, more stars can form, raising the average SFR in all galaxies. This could be counterbalanced by a lower SFR per gas mass (SFR/Mgas) at higher redshift but existing works have shown that SFR/Mgas is either relatively constant or rises with redshift (Tacconi et al. 2010; Scoville et al. 2016).

It has been argued that the slope of the MS should be unity for all mass ranges once the SFR has become stable in a galaxy (Pan et al. 2017). However, the slope of the MS is typically found to be lower than unity, with values between 0.4 and 1.0 (e.g. Whitaker et al. 2012; Speagle et al. 2014; Tomczak et al. 2016). This discrepancy is believed to be the result of a combination of quenching and a reduction in the relative size of the cold gas reservoir as the M⋆ of a galaxy increases (Pan et al. 2017).

Far-infrared (FIR) emission is a key component for accurately determining the SFR of an object. Part of the ultraviolet (UV) emission from young stars heats dust in a galaxy, which then re-radiates in the infrared (IR). Observing in the UV would provide a direct measure of SFR, but could underestimate the total SFR as this absorption by dust obscures some of the UV emission (e.g. Meurer et al. 1999; Dale et al. 2009; Bourne et al. 2017; Dunlop et al. 2017). As a result, observing IR emission in combination with UV emission is key to providing a complete picture of the SFR.

However, FIR surveys, such as those conducted with the ESA Herschel Space Observatory (Pilbratt et al. 2010) Spectral and Photometric Imaging Receiver (SPIRE; Griffin et al. 2010), have relatively poor resolution with respect to optical, and also have source confusion (e.g. Nguyen et al. 2010; Oliver et al. 2012; Hurley et al. 2017). Previously, stacking was used to recover the flux densities of fainter sources (e.g. Pannella et al. 2015; Schreiber et al. 2015; Harris et al. 2016; Álvarez-Márquez et al. 2016). However, by its very nature, stacking only provides the average properties (e.g. mean, median) of the objects in the stack1. As a result, the properties of individual objects cannot be determined, which results in the loss of information about the wider properties of a galaxy population.

To overcome confusion, it is necessary to de-blend the maps to generate individual flux density measurements for both faint and bright sources. For SPIRE, this was done with the De-blended SPIRE Photometry algorithm (DESPHOT; Roseboom et al. 2010, 2012; Wang et al. 2014), amongst other techniques (e.g. Laidler et al. 2007; Béthermin et al. 2010b; Kurczynski & Gawiser 2010; Viero et al. 2013; Merlin et al. 2015; Safarzadeh et al. 2015; Wright et al. 2016). DESPHOT misassigns flux densities when more than one source is within a beam, and is also unable to realistically derive flux density errors of a given source (see discussion in Hurley et al. 2017). To improve source de-blending, XID+ (Hurley et al. 2017, see also Sect. 3.1.2) has been developed and subsequently expanded to improve flux density estimation by including more precise flux density priors (Pearson et al. 2017).

With XID+, it has been shown that more than 95% of blindly detected SPIRE 250 μm sources witha flux density greater than 30 mJy contain more than one object which contribute more than 10% of the source’s total flux density. At least 70% of the flux density from these sources is assigned to the two brightest objects (Scudder et al. 2016). This suggests that many current SPIRE catalogues have too much flux density assigned to their objects, and any derived physical parameter that relies on the SPIRE emission, such as the total infrared luminosity or star formation rate (SFR), will be overestimated.

In this work, XID+ is used to de-blend the maps in the SPIRE bands, resulting in a catalogue of over 200 000 objects with SPIRE flux density measurements in the COSMOS field, compared to approximately 32 000 sources blindly detected at 250 μm (Roseboom et al. 2010; Wang et al. 2014; Sparre et al. 2015). These data, along with UV to Herschel Photodetector Array Camera andSpectrometer (PACS; Poglitsch et al. 2010) data from a multi-wavelength catalogue, are then used to examine the MS.

The paper is structured as follows. Section 2 describes where the data were collected. Section 3 explains the methodology and tools used to de-blend the SPIRE bands and find the MS. Section4 explores the results of the analysis. Section 5 provides discussion, while Sect. 6 provides a summary of the conclusions. Where necessary, Wilkinson Microwave Anisotropy Probe year 7 (WMAP7) cosmology (Komatsu et al. 2011; Larson et al. 2011) is followed, with ΩM = 0.272, ΩΛ = 0.728, and H0 = 70.4 km s−1 Mpc−1, to be consistent with CIGALE2 (Noll et al. 2009, see also Sect. 3.1.1).

2 Data sets

For this work, the COSMOS field (Scoville et al. 2007) was chosen due to the wealth of multi-wavelength data available within its two square degree coverage. It also benefits from Herschel SPIRE data, which is considered to be homogenous and of high quality (e.g. Pearson et al. 2017).

For the multi-wavelength data, the public COSMOS2015 catalogue (Laigle et al. 2016) was used with CIGALE to generate flux density priors in the 250 μm, 350 μm, and 500 μm SPIRE bands for XID+ (see Sects. 3.1.1 and 3.2). COSMOS2015 contains photometric data in over 30 bands, narrow, medium, and broad, for approximately 1.2 × 106 objects in the COSMOS field. Here, we only use the bands that cover the UV to PACS 160 μm in CIGALE; a full list of the bands used can be found in Table 1. The aperture photometry was converted to total photometry and corrected for foreground extinction following Laigle et al. (2016). As we wanted to use the Multi-band Imaging Photometer for Spitzer 24 μm (MIPS24; Rieke et al. 2004) data for the ≈40 000 objects with MIPS24 data to help constrain the FIR SED, we used objects that only fall within the MIPS24 image in COSMOS. The MISP24 coverage falls within the SPIRE coverage; the requirements are as follows (see also Fig. 1): (2)

The latest SPIRE images from the Herschel Database in Marseille3, the Data Release 4 maps from the Herschel Multi-tiered Extragalactic Survey (Oliver et al. 2012), were used in XID+ to extract the SPIRE flux densities. The 250 μm, 350 μm, and 500 μm band images have beam full widths at half maximum of 18.1′′, 25.2′′, and 36.6′′ (Griffin et al. 2010)and 5σ confusion limits of 24.0, 27.5, and 30.5 mJy, respectively (Nguyen et al. 2010).

Table 1

Telescopes and associated bands that were used for CIGALE spectral energy distribution fitting from the COSMOS2015 catalogue.

thumbnail Fig. 1

Imageof the SPIRE 250 μm COSMOS coverage (Oliver et al. 2012, 7.84 deg2, red) with the overlayed MIPS 24 μm COSMOS coverage (Sanders et al. 2007, 2.23 deg2, blue). The data used in this work were cut to match the MIPS 24 μm coverage.

3 Methodology

3.1 Tools

3.1.1 CIGALE

Code Investigating GALaxy Emission (CIGALE; Noll et al. 2009) is a spectral energy distribution (SED) modelling and fitting tool with an improved fitting procedure by Serra et al. (2011). Here, the Python version 0.11.0 is used (Burgarella et al. 2005; Noll et al. 2009, Boquien et al., in prep.) to generate SEDs and to fit these SEDs to the UV to PACS data from COSMOS2015 to estimate the SPIRE 250 μm, 350 μm, and 500 μm flux densities for use as a flux density prior in XID+ (see Sects. 3.1.2 and 3.2). After the SPIRE band flux densities had been extracted, CIGALE was also used to calculate the physical parameters of each object, such as SFR and M⋆ as well as rest frame colours. CIGALE uses the energy balance between the attenuated UV emission by dust and the IR emission, allowing the estimation of the FIR flux densities. The reported values and errors for the SPIRE flux densities and physical parameters are created using CIGALE’s Bayesian probability density function analysis. CIGALE also gives the flux densities and physical parameters of the best fitting model for each object, but these are not used in this work.

The choices for the SED model components and parameters for the SPIRE band priors follow Pearson et al. (2017), with a different dust attenuation model and other minor changes, and will briefly be repeated here. We use a delayed exponentially declining star formation history (SFH) with an exponentially declining burst, Bruzual & Charlot (2003) stellar emission, Chabrier (2003) initial mass function (IMF), Charlot & Fall (2000) dust attenuation, the updated Draine et al. (2014) version of the Draine & Li (2007) IR dust emission, and Fritz et al. (2006) AGN models. The dust attenuation model was changed from Calzetti et al. (2000), used in Pearson et al. (2017), to Charlot & Fall (2000) as recent work by Lo Faro et al. (2017) has shown that the Calzetti et al. (2000) dust model cannot accurately reproduce the attenuation seen in a sample of dusty galaxies at z ≈ 2. A list of parameters used, where they differ from the default values, along with a justification can be found in Appendix A.

3.1.2 XID+

XID+4 (Hurley et al. 2017) is a probabilistic de-blending tool used to extract source flux densities from photometry maps that suffer from source confusion. This is achieved by using Bayesian inference to explore the posterior. Once converged, the flux density is reported along with the upper and lower 1σ uncertainties. In its standard form, XID+ uses a flat prior in parameter space, between zero and the brightest value in the map, along with the source positions on the sky.

This work follows Pearson et al. (2017) by using a more informed Gaussian prior, again truncated between zero and the brightest value in the map. The mean and sigma for these Gaussian priors are generated by using CIGALE models to estimate the flux densities for the mean and using two times the error on these estimates as the sigma. To allow parallelisation, which reduces the time taken for XID+ to de-blend the map, the map is split up into tiles based on the Hierarchical Equal Area isoLatitude Pixelization of a sphere system (HEALPix; Górski et al. 2005) using order 11, which corresponds to an area of2.95 arcmin2 per tile. Order 11 was chosen as it is a compromise between the number of objects in a tile (more objects means a more reliable flux density extraction) and the time it takes a tile to run.

3.2 Science pipeline

The extraction of the flux densities in the SPIRE bands, which follows a similar pipeline to the one used in Pearson et al. (2017), is briefly summarised in Fig. 2 and is repeated here for completeness. We begin by using the far-UV (FUV) to PACS 160 μm data from COSMOS2015 to generate estimates for the flux densities in all three SPIRE bands simultaneously using CIGALE’s Bayesian analysis. All objects that were not classified as galaxies (TYPE flag in COSMOS2015 not set to 0, approximately 1.3%) were removed, as were objects without any photometric redshift (ZPDF ≤ 0, approximately 2.9%). All predicted flux densities were then used in XID+ to extract the flux densities for the objects; we did not remove any objects with poor χ2 values. Once the SPIRE flux densities were extracted, these SPIRE data were added to the FUV-PACS data and CIGALE was rerun, this time to get estimates for M⋆ and SFR. The same CIGALE models were used for the flux estimation and to obtain the physical parameters so the results from each CIGALE run should not be degenerate.

The first CIGALE run provides flux density estimates for all the objects in the catalogue. However, the flux density estimates for the faintest objects will be highly uncertain and hence result in extracted flux densities that are unreliable. To find the depth to which we can reliably run XID+, a number of different cut depths on the predicted flux densities at 250 μm were used, from 0.2 mJy to 20 mJy, and XID+ run on 25 tiles. A residual map was created by subtracting the replicated map from XID+ from the original image for each depth. In the ideal situation where all sources are correctly accounted for, the residual map should have a scatter consistent with 1σ instrument noise. Figure 3 shows how the scatter of the residual 250 μm map changes with cut depth with and without 3σ clipping. As can be seen, as the cut gets deeper, from 20 mJy to 1 mJy, the scatter reduces; cut depths deeper than 1 mJydo not cause much of a change. All cut depths have a scatter greater than the instrument noise at 250 μm (1.71 mJy). This is to be expected as our source list will miss the faintest sources and the flux densities assigned to the known sources will not exactly coincide with the true flux densities for all sources. As a result of all the deeper cuts having a scatter larger than the 1σ instrument noise, and wanting to keep as large a sample as possible while not risking creating too many degeneracies, a predicted 250 μm cut of 0.7 mJy was used. For this cut, we ignore any uncertainty on the flux density estimate. This leaves 205 958 objects to be run through XID+.

To prevent an overly restrictive prior in XID+ and conservatively capture the uncertainty in SED modelling, the errors on the flux density estimates are expanded by a factor of two (see also Appendix B). The estimates of all three SPIRE bands from CIGALE are then used as the priors for XID+. The flux density estimates are used as the means in the XID+ priors while the expanded errors are used as the standard deviations. XID+ was then run on the SPIRE images. The next step is to run the data from COSMOS2015 and XID+ through CIGALE to generate the required parameters for each object: M⋆, SFR, and U–V and V–J colours; the colours are needed for the QG/SFG separation.

thumbnail Fig. 2

Brief summary of the science pipeline.
thumbnail Fig. 3

Scatter of the 250 μm residual map using different depth cuts on the prior list. The blue and red points are with and without 3σ clipping, respectively, while the green line is the 1σ instrument noise for the COSMOS field.

3.3 Masscompleteness and redshift binning

To study the MS, we removed all objects with a redshift below 0.2 as the mean error on the photometric redshifts for objects below z = 0.2 is approximately 0.2 (for discussion on the redshifts and their errors, see Laigle et al. 2016). The remaining objects were then binned by redshift. Each bin’s width was determined by using the average error on the redshift of the objects within a bin; the bin size was the average error of all the objects in the bin rounded up to the next decimal place. Using a bin size of twice the error was also briefly explored and provided consistent results.

Once binned by redshift, the galaxies were cut for completeness. This completeness cut was done empirically by following Pozzetti et al. (2010). The Ks band magnitude limits (Ks lim) were set to 24.7, the 3σ limit for the objects within the UltraVISTA (McCracken et al. 2012) ultra-deep stripes, and 24.0, the 3σ limit for the rest of the COSMOS field, in the ultra-deep and deep regions of COSMOS, respectively (Laigle et al. 2016). For each galaxy with a redshift (ZPDF > 0) and a Ks band detection the mass the galaxy would need (Mlim) to be observed at the magnitude limit was calculated using (3)

where M is the galaxy’s mass in M⊙ and Ks is the galaxy’s Ks band magnitude. In each redshift bin, the faintest 20% of objects were selected and the limiting mass was the Mlim value which 90% of these faintest objects lie below. As not all objects have a Ks band detection (the catalogue is zYJHKs selected), these limiting mass data points were used for all objects. The mass limits for the deep and ultra-deep regions of COSMOS are shown in Fig 4.

thumbnail Fig. 4

Masses of all objects detected in the Ks band (blue) are shown against redshift along with the faintest 20% in each redshift bin (red) for the (panel a) deep and (panel b) ultra-deep regions. The 90% completeness limit is shown by the thick black lines, while the dashed black lines show the edges of each redshift bin.

3.4 Forward modelling

To determine the MS in each redshift bin, we use the Markov chain Monte Carlo (MCMC) sampler emcee (Foreman-Mackey et al. 2013)to sample the parameter space of the chosen MS model. For our routine, we create model SFRs using the observed M⋆, observed redshift, and the MS being tested at each step. For each M⋆, a random number is drawn from a Gaussian distribution centred on the SFR of the MS at that M⋆ and with the scatter of the MS as the standard deviation. The Gaussian distribution is also truncated such that it reproduces the observed upper and lower SFR limit. These upper and lower SFR limits are determined by finding the SFR in each redshift bin that 0.1% of the objects fall above or below, and then fitting to these values as a function of redshift using (4)

where S is log (SFR∕M⊙ yr−1) and Bn are the coefficients that are found. The parameters B0, B1, and B2, take the values 0.61 ± 0.32, –0.37 ± 0.01 M⊙ yr−1, and 3.13 ± 0.13 log(M⊙ yr−1) for the upper limit and 2.86 ± 0.19, –0.17 ± 0.05 M⊙ yr−1, and 0.12 ± 0.10 log(M⊙ yr−1) for the lower limit, respectively. The upper and lower limits are applied to each simulated object individually using the observed redshifts.

Once the model SFRs have been generated, both the SFR and M⋆ are perturbed by adding a second random number drawn from a Gaussian centred on zero and with a standard deviation equal to the error on the observed SFR or M⋆ of that individual object. An example of the MS generated at one step, assuming a linear power law, is shown in Fig. 5, along with the MS used to generate at that step and the contours of the observed data.

To compare the model data to the observed data at each step, the two data sets are binned by stellar mass into identical bins with a width of 0.25 dex. The mean and standard deviation of the SFRs in each mass bin are calculated and the mean and standard deviations of the model data are compared to their counterparts from the observed data. The smaller the differences, the greater the likelihood that the model is a correct representation of the observed data.

All the parameters in the models used were treated as if they were uncorrelated, although this is not strictly true. However, this method was found to accurately recover input relations used to generate mock sets of data that assume a linear power law. Figure 6 shows an example of the posterior for a mock data set with known slope (0.6), normalisation (0.7 log (M⊙∕yr−1)), and scatter (0.3 dex). As can be seen, the input parameters are recovered within error.

thumbnail Fig. 5

Example of the data generated at one step of the MCMC routine for the lowest (0.2 ≤ z < 0.5) redshift bin, shown as number density from high (dark red) to low (dark blue). The MS being tested at this step, in this case the most likely step, is shown as the red line, while the contours for number density of the observed data are shown in orange. The size of the average observed error on SFR and M⋆ is also shown as a blue cross.
thumbnail Fig. 6

Corner plot (Foreman-Mackey 2016) of the marginalised posterior of the forward modellingroutine applied to a mock data set with known slope (0.6), normalisation (0.7 log (M⊙∕yr−1)), and scatter (0.3 dex), shown as the blue lines. Panels on the diagonal show the 1D marginalised posteriors for the slope, normalisation and scatter (left to right panels). Off-diagonal panels show the combined 2D posteriors as labelled by their axes. The recovered 16th, 50th, and 84th percentiles are shown by the dashed vertical lines; all the input parameters are recovered, within error.

4 Results

The M⋆ found using CIGALE are compared to those found by Laigle et al. (2016) in COSMOS2015. We find that our M⋆ are higher, on average, by 0.15 dex and are consistent with COSMOS2015 within the COSMOS2015 average error of 0.15 dex, but not the 0.10 dex average error of our M⋆. Both Laigle et al. (2016) and this work use the Chabrier (2003) IMF and the Bruzual & Charlot (2003) stellar population model, so this is not the cause of the slight discrepancy, but the dust attenuation models differ: this work uses the Charlot & Fall (2000) dust attenuation, while Laigle et al. (2016) uses the Calzetti et al. (2000) dust attenuation. Recently, Lo Faro et al. (2017) have investigated the effects of dust attenuation laws on the physical properties of dusty galaxies at z ≈ 2. They find that using the Charlot & Fall (2000) dust attenuation results in higher M⋆ values than Calzetti et al. (2000), resulting in M⋆ higher by a factor of 1.4, or approximately 0.15 dex. This is consistent with the difference in M⋆ found between our current work and Laigle et al. (2016). Thus, we conclude that the difference in M⋆ is a result of thechoice of dust attenuation model.

To remove the QG, we follow the Whitaker et al. (2011) UVJ colour cut (5)

where the rest-frame (U–V) and (V–J) colours come from the second CIGALE run. If these conditions are not met, the object is SF. The (U–V) and (V–J) criteria for 2.0 < z < 3.5 were expanded for all objects with a redshift greater than 2, such that the final line in Eq. (5) becomes (6)

Once the QG objects are removed, we fit two models to the data: the Lee et al. (2015) description, which contains a turn-over (7)

where S is log (SFR∕M⊙ yr−1), S0 is the maximum value of S that the function approaches at high mass, M0 is the turn-over mass in M⊙, and γ is the low-mass slope, and the Whitaker et al. (2012) single power law (8)

where α is the slope and β is the normalisation at log(M⋆∕M⊙) = 10.5.

Fitting with Eq. (7) was reasonably unsuccessful. No redshift bins show significant evidence of a high-mass turn-over with all redshift bins being consistent with a simple power law, as can be seen in the example in Fig. 7. The majority of the turn-over masses found from the MCMC forward modelling are larger than the highest mass object within each redshift bin, approximately log (M⋆∕M⊙) = 12.5. As a result, these cannot be considered reliable results as the turn-over position is unconstrained by the data. The turn-over positions at all redshifts, assuming they are reliable, are also much higher than those found in the literature (log (M0∕M⊙) ≈ 10.5 e.g. Lee et al. 2015; Tomczak et al. 2016). We also found that when fitting Eq. (7) with γ held fixed at 1.17, which was shown by Lee et al. (2015) to be reasonable and also leaves the same number of free parameters as Eq. (8), the majority of the redshift bins have a best fit that is less probable than the best fit using Eq. (8). As a result, we conclude that Eq. (8) is a better description of our data.

With the knowledge that our data do not support a turn-over in the MS, Eq. (8) was also fitted to our data, as shown in Fig. 8. This form of the MS was much more successful, with well-constrained fitting parameters in all redshift bins (see Table 2 and Fig. 9). The β parameter (normalisation) increases rapidly with redshift out to z ≈ 2 before increasing more slowly (see magenta points in Fig. 10b). The evolution of β was found, fitting using the SciPy (Jones et al. 2001) orthogonal distance regression (ODR; Boggs & Rogers 1990) package, to follow (9)

with errors in the coefficients from the ODR fitting.

The α parameter (slope) is slightly less smooth in its evolution (see Fig. 10a, magenta). Generally, the slope increases with redshift across the entire redshift range of this study. There is a potential rise between z ≈ 1.8 and z ≈ 2.9, although this is very marginal and may be a result of the SED model gridding causing a higher slope (see Fig. 8). A linear relation was used to find the slope evolution, which was determined to be (10)

Short tests were also done using just the deep and ultra-deep data. At redshifts below 4.9, just using the deep or ultra-deep data made no significant difference to the α and β parameters: the values change within error, compared to using both the deep and ultra-deep data together. With just the deep data, the slope of the highest redshift bin becomes much lower, although it is consistent within error. For normalisation, the deep data has a lower normalisation, but again it is consistent within error. This is likely a result of low number statistics. There are only 16 UVJ selected objects with z ≥ 4.9 in the deep coverage. With so few objects to fit with, the results are highly uncertain with just the deep data.

After determining the MS, we check the scatter to ensure that we have actually found the MS. The intrinsic scatter was determined during the MCMC fitting and was assumed to be independent of M⋆, which is believed to be true (Whitaker et al. 2012; Speagle et al. 2014; Tomczak et al. 2016). Our intrinsic scatter is found to be consistent with other literature (e.g. Whitaker et al. 2012; Speagle et al. 2014; Sparre et al. 2015; Schreiber et al. 2015; Tomczak et al. 2016, see Fig. 11), suggesting that we have indeed found the MS. We also fit a flat line to our intrinsic scatter, assuming that there is no redshift evolution, which gives an amplitude of 0.23± 0.03 dex, lower than the expected 0.3 dex typically found in observational works but consistent with the results from the Illustris Simulation (Genel et al. 2014; Vogelsberger et al. 2014; Sparre et al. 2015). Our small intrinsic scatter may be an underestimation resulting from the M⋆ and SFR being derived from the same SED fitting procedure. This results in some correlation between these two quantities and hence reduces the scatter about the MS.

thumbnail Fig. 7

Comparison of fitting with Eq. (7) (orange dashed line) and Eq. (8) (red line) in the 0.5 ≤ z < 0.8 redshift bin. The galaxies are shown as a number density plot, with dark red being high number density and dark blue low number density, and the size of the average error on SFR and M⋆ is shown as a blue cross. There is very little difference in shape between the two fits, demonstrating that Eq. (8) is the preferred form of the MS.
thumbnail Fig. 8

Fitting of Eq. (8) to the objects in the redshift bins as labelled. The solid line is the most likely MS across the fitted M⋆ range, while the dotted line is an extrapolation across the M⋆ range of theplot. The galaxies are shown as a number density plot, with dark red being high number density and dark blue low number density. The vertical density discontinuities are the result of the two depths of data used: the deep mass limit is the yellow dashed line and the ultra-deep mass limit is the dot-dashed yellow line. Each panel also indicates the χ2 of the most likely MS () and the number of degrees of freedom in the fitting (ndof), and shows the size of the average SFR and M⋆ errors as a blue cross.
Table 2

Parameters from fitting Eq. (8) to the UVJ selected star forming galaxies, with α the slope and β the normalisation at log(M⋆/M⊙) = 10.5.

thumbnail Fig. 9

Fitted MS trends using Eq. (8). The solid and dashed lines are the MS across the fittedM⋆ range. The dotted lines are extrapolations across the M⋆ range of the plot. The normalisation of the MS clearly increases as redshift increases. The slight decrease in slope at low redshift can be seen along with the increase above z = 1.1. The very steep slopes at 1.8 ≤ z < 2.9 can also be seen.
thumbnail Fig. 10

Comparison of the UVJ selected MS results of this work (magenta) with the observational MS from Speagle et al. (2014, FUV data in blue; IR data in green; and radio, hydrogen lines, and UV SED fitting in yellow), Schreiber et al. (2015) low mass MS (red), and the MS from the Illustris Simulation (Vogelsberger et al. 2014; Sparre et al. 2015, black diamonds). The α and β parameters from Eq. (8) are in panels a and b, respectively. As Schreiber et al. (2015) hold their low-mass slope constant at unity, this is indicated in panel a as a flat red line. The redshifts shown for this work are the mean redshift in each redshift bin, while the horizontal error bars show the width of the redshift bin. A version of this plot using SFRs derived from IR template fitting can be found in Appendix C.
thumbnail Fig. 11

Intrinsic scatter about the MS found from the MCMC fitting of the data in each redshift bin to Eq. (8). This work (magenta) is shown along with the observed scatters from Speagle et al. (2014, FUV data in blue, IR data in green, and other data in yellow) and Schreiber et al. (2015, red) as well as the median scatter found at each redshift in the Illustris Simulation (Genel et al. 2014; Vogelsberger et al. 2014; Sparre et al. 2015, black diamonds). The redshifts shown for this work are the mean redshift in each redshift bin while the horizontal error bars show the width of the redshift bin. Also shown is the best fit to this work’s intrinsic scatter assuming no redshift evolution (dashed magenta line). The intrinsic scatter found in this work is consistent with existing literature, although above z ≈ 1.8 our intrinsicscatter is smaller than is found in works that use IR SFR traces.

5 Discussion

5.1 Formof the main sequence

In this work, we find little evidence of a turn-over. Thus, any cause of the reduction in specific SFR (SFR/M⋆, sSFR) with mass does not suddenly ‘turn on’, rather it is an effect that slowly increases with M⋆. An increasing slope with redshift indicates that the mass dependent reduction in sSFR becomes weaker as we look further back, meaning that the mechanism behind this was less strong in the past than it is now.

Second, at relatively low redshift we find that the slope is relatively shallow, which suggests that the ability of galaxies to form stars reduces as mass increases. This may be a result of a decrease in star formation efficiency, a reduction in the available cold gas, or some other mass-dependent quenching mechanism. We also see that the slope increases with redshift, implying that galaxies in the younger universe were more self-similar than they are today. Thus, in the modern universe, large galaxies have lower specific SFRs (sSFR) than smaller galaxies, while at high redshift, the sSFRs of large and small galaxies are similar.

Third, we find that the slope of the MS is very high between redshifts of 1.8 and 2.9 with respect to other epochs. We do note, however, that the significance of this rise is very marginal. A steep slope at these times suggests a higher sSFR in large galaxies compared to other times. This increase in the slope coincides with the peak of the cosmic star formation history (CSFH) at z ≈ 2−3 (Madau & Dickinson 2014). It may also simply be a result of the SED model gridding causing a higher slope, indicated by the discontinuities in the two top right panels of Fig. 8.

When the normalisation is also considered, the low-mass end of the MS also sees changes in a galaxy’s ability to form stars. Eq. (8) has a normalisation at log(M⋆∕M⊙) = 10.5 and it does not evolve fast enough to maintain a constant crossing point of the MS at all redshifts; the intercept mass of the MS in one redshift bin with the MS in the lowest redshift bin increases out to z ≈ 1.8 before becoming approximately constant, as can be seen in Fig. 9. Thus, the low-mass galaxies must have a decrease in sSFR as redshift increases, out to z ≈ 1.8, while the high-mass galaxies have an increasing sSFR.

5.2 Contamination by quiescent galaxies

We find the main sequences in the highest redshift bins (z > 1.8) are identical regardless of whether we cut for SFG or not. This may be a result of the dusty SFG and QG at high redshift having their near-IR (NIR) observed frame emission suppressed below the detection limits of the NIR selected catalogue by dust obscuration in the galaxy. Alternatively, there is a greater mixing of SFG and QG at lower redshifts and higher masses, pulling down the MS when the QG classed objects are removed. The latter is likely the correct scenario. There is little evidence for a turn-over at lower redshifts after selection of the SFGs. We thus surmise that the turn-over in other works may be the result of a SFG selection that is not effective enough in removing the QG. This would result in an overabundance of low SFR galaxies at higher mass which act to pull the high-mass MS down, similar to the conclusions of Johnston et al. (2015). It is the lack of turn-off in the MS at lower redshifts that leads us to believe that the high-z colour cut employed (Eq. (5)) can be extended to all redshifts above z = 2.0 as we also recover the linear MS at these higher redshifts.

5.3 Comparisons with other works

The work by Speagle et al. (2014) collates a large number of MS studies and generates a consensus regarding the MS. This provides an excellent base line to which we can compare this work. As with the consensus result, this work finds little evidence of a high-mass turn-over. Broadly, our results agree with Speagle et al. (2014): both the slope and normalisation increase with redshift; the rate of change in normalisation decreases as redshift increases. However, quantitatively there are differences.

The normalisation found in this work is in reasonable agreement with the Speagle et al. (2014) collection (see Fig. 10b), but lies towards the bottom of the collection. This slight lowering of the normalisation, especially compared with studies using IR SFR indicators, may be due to the splitting of brighter objects into a number of fainter sources.

The slope of this work is consistent with the Speagle et al. (2014) results at all redshifts, rising from the bottom of the collection at low redshifts to the top of the collection at high redshifts. Between z ≈ 2 and z ≈ 3 we are approximately consistent with the far-ultraviolet SFR indicated MS, but are typically higher than the MS using other (radio, hydrogen lines, and UV SED fitting) SFR indicators. This is potentially a result of the data in this work looking through the dust with the Herschel data, which will raise the SFR of the high-mass, dusty galaxies. However, this should also result in the FUV MS having a shallower slope than this work as well.

The work by Schreiber et al. (2015) looked at the MS using Herschel, but their study relied on stacking rather than de-blending.We convert the Schreiber et al. (2015) results from the Salpeter (1955) IMF to the Chabrier (2003) IMF used here, following Speagle et al. (2014, see Eq. (2)) to correct the mass and Kennicutt & Evans (2012) to correct the SFR. Their results found a turn-over at high mass, unlike the results of this study. This high-mass turn-over is not a result of stacking; we performed our own stacking analysis and found no turn-over or UVJ selection because our UVJ selection is similar to that of Schreiber et al. (2015; see Appendix C).

As a result of their turn-over, just the low-mass slope of Schreiber et al. (2015) is used for comparison in this work. Schreiber et al. (2015) also enforce a slope of unity while the normalisation is allowed to change. A non-redshift dependant slope of unity is not found to be a good description of our data, with the slope remaining below this until the highest redshift bin. We also find that our normalisation is considerably lower at all redshifts, by approximately 0.4 dex, with the difference being larger at lower redshifts (see Fig. 10). However, this difference in normalisation is a result of the different methods used to derive the SFR of the galaxies, with our SED-based SFR producing values lower by 0.4 dex when compared to the SFRs derived from the IR luminosity from CIGALE. We also see the same 0.4 dex offset when deriving our own IR luminosity SFRs by fitting our SPIRE data to Chary & Elbaz (2001) templates (see Appendix C).

We also compare our results to those from the Illustris Simulation (Genel et al. 2014; Vogelsberger et al. 2014; Sparre et al. 2015). We fit Eq. (8) to the SFR-M⋆ positions in Fig. 1 of Sparre et al. (2015) for the z = 0, 1, 2, and 4 redshift bins. The z = 0 bin shows clear evidence of a turn-over, but we still fit a simple power law to all the data in this redshift bin. The resulting α and β values can be found in Fig. 10 as black diamonds. Our work has slopes that are much shallower than the Illustris values below a redshift of ≈ 1.8. At z ≈ 2, this work and Illustris agree withinerror, while at z ≈ 4 our slope appears to be consistent with Illustris. The normalisation of the MS from the Illustris simulation, along with other simulations, is known to be too low with respect to observations, certainly at z = 1−2 (Sparre et al. 2015; Furlong et al. 2015). Our normalisation results are not consistent with this picture; the normalisation we find is consistent with the results from Illustris. This may suggest that the low normalisations found in simulations may be the result of simulations not being affected by the source blending that observational studies can show.

6 Conclusions

The MS is a tight relation between M⋆ and the SFR of galaxies in the universe. The FIR is a very important part of the spectrum for determining the SFR of an object; however, current cosmological FIR surveys suffer from poor resolution. Here, the latest de-blending tool, XID+, and techniques have been used to break through the Herschel confusion limit, allowing better SFR estimates to be generated for nearly a quarter of a million objects in the COSMOS field. This has allowed the examination of the main sequence of star forming galaxies beyond the Herschel confusion limit.

The objects were binned by redshift between 0.2 ≤ z < 6.0 and fitted with both a power law and a power law with a high-mass turn-over. The high-mass turn-over form of the MS failed to provide reasonable values for its parameters. The simple power law, on the other hand, produced well-behaved parameters in every redshift bin. Thus, we conclude that the MS using multi-wavelength data and our de-blended Herschel SPIRE data shows little evidence for a high-mass turn-over.

We find that the normalisation of the MS increases with redshift, rapidly out to z ≈ 2 and more slowly thereafter. The slope has a different evolution, rising steadily across all redshifts, with perhaps a weak indication of a brief peak at 1.8 ≤ z < 2.9. An increase in slope with redshift is likely a result of high-mass galaxies having an increase in star formation efficiency with respect to the low-mass galaxies. This results in all galaxies becoming more self similar with increasing redshift, as the slope increases towards unity. The brief peak in slope, if present, appears to coincide with the peak of the CSFH.

Our normalisation results are consistent with the Speagle et al. (2014) MS compilation results within error, and with the Schreiber et al. (2015) stacked Herschel MS, once adjusted for IMF and SFR tracers. The slope of our MS is also consistent with the Speagle et al. (2014) compilation but moves from the lower regions at low redshift to the upper regions at higher redshifts. This may be a result of hidden star formation being revealed by the Herschel data, demonstrating the importance of the FIR/sub-mm for the MS.

Acknowledgements

We would like to thank the anonymous referee whose comments greatly improved this paper. We would like to thank the Center for Information Technology of the University of Groningen for their support andfor providing access to the Peregrine high performance computing cluster. This research has made use of data from HerMES project. HerMES is a Herschel Key Programme utilising Guaranteed Time from the SPIRE instrument team, ESAC scientists, and a mission scientist. The HerMES data were accessed through the HerschelDatabase in Marseille (HeDaM - http://hedam.lam.fr) operated by CeSAM and hosted by the Laboratoire d’Astrophysique de Marseille. HerMES DR4 was made possible through support of the Herschel Extragalactic Legacy Project, HELP. PDH and SJO would like to acknowledge the research leading to these results, which has received funding from the European Union Seventh Framework Programme FP7/2007-2013/ under grant agreement n° 607254. This publication reflects only the author’s view and the European Union is not responsible for any use that may be made of the information contained therein.

Appendix A CIGALE Parameters

In this workwe closely follow the CIGALE model selection of Pearson et al. (2017, Appendix A). A delayed exponentially declining SFH with an exponentially declining burst.

was used over the more commonly used exponentially declining or delayed exponential SFH as these two SFHs did not appear to reproduce the expected starburst population in the star formation rate versus stellar mass plane. The e-folding time of the two stellar populations (old and young) in the SFH was roughly matched to that of Mitchell et al. (2013). As Mitchell et al. (2013) used a single declining exponential SFH, the e-folding times were split with the burst population taking values of 9 Gyr and above, and the main population taking values of less than 9 Gyr. For the ages of the main population, the values were sampled linearlybetween 1 and 11 Gyr with a wider sampling up to 13 Gyr. The mass fraction of the burst population follows Ciesla et al. (2015) along with the age of the young stellar population, which also had a lower age of 0.001 Gyr added. The Bruzual & Charlot (2003) stellar population model was used with a Chabrier (2003) initial mass function. Asthis study was not to explore the metallicity of galaxies, it was decided to leave the metallically at solar.

Recent work by Lo Faro et al. (2017) has shown that a single power law is a poor representation of the attenuation that occurs in dusty galaxies. As such, we adopt a double power law, similar to that of Charlot & Fall (2000), with individual power laws for the birth clouds (BCs) and inter stellar medium (ISM). However, we use a smaller range of values for the V-band attenuation in the BCs than Charlot & Fall (2000), a maximum value of 3.8, as tests with higher attenuation values do not noticeably effect the results.

For the dust emission, the polycyclic aromatic hydrocarbon (PAH) fraction had an increase in range around the default 2.5 so more fractions could be sampled while keeping the number of models reasonable. The minimum scaling factor of the radiation field was similarly given a range to sample with an increase in the smallest value from 1.0 to 5.0. The illuminated fraction was reduced to 0.02, following Ciesla et al. (2015).

The parameters in the AGN model were matched to those used by Ciesla et al. (2015), who undertook adetailed study of AGN host galaxy emission using CIGALE. The number of choices of fracAGN was reducedfrom 14 to 5, while still covering the same range, to reduce the number of models created by CIGALE and hence decrease runtime.

A list ofparameters, where they differ from default, can be found in Table A.1.

Table A.1

Parameters used for the various properties in the CIGALE model SEDs where they differ from the default values.

Appendix B XID+ prior width

We require the prior from CIGALE to be as informative as possible in XID+ while not being overly restrictive. To determine the preferred expansion factor of the SPIRE flux density estimate errors from CIGALE, we ran XID+ on the tiles that contain at least one of 178 objects detected at 870μm from Scoville et al. (2014, 2016) three times; expanding the flux density estimate error by a factor of 2, 3, and 4. The SPIRE flux densities from XID+ were then added to the data from COSMOS2015 and CIGALE was rerun to predict the Atacama Large Millimeter/submillimeter Array (ALMA) 870 μm flux densities, for each of the three XID+ results. These predictions were compared to the ALMA 870 μm observations by subtracting the CIGALE predictions from the observations and fitting a Gaussian to the resulting distribution. The histograms of the comparisons, along with the fitted Gaussian distributions, are found in Fig. B.1. As the standard deviations are all approximately consistent, we chose to expand the errors from CIGALE by a factor of two as this provides the smallest deviation of the mean from zero.

thumbnail Fig. B.1

Distributions for the difference between the observed ALMA 870 μm flux densities and the predicted 870 μm flux densities from CIGALE using an error expansion factor of 2 (panel a), 3 (panel b), and 4 (panel c) in XID+. The Gaussian distribution for each expansion factor is also shown in orange, along with the means (μ) and standard deviations (σ) of the distributions. All the distributions have an approximately consistent σ so the best expansion factor was deemed to be that with μ closest to zero: 2 times the error from CIGALE.

Appendix C Further comparisons with Schreiber et al. (2015)

As mentioned in Sect. 5.3, the normalisation presented in this work is consistent with that of Schreiber et al. (2015) once we account for IMF and SFR tracers. We also note that the high-mass turn-over found in Schreiber et al. (2015) is not a result of stacking. We discuss this here in further detail.

The disparity in normalisation is likely a result of the different methods of determining SFR. Schreiber et al. (2015) convert the UV and IR luminosities into UV and IR SFRs and then combine these values to generate a total SFR, while we use full SED fitting. As a result, our SFRs are lower as not all of the IR luminosity will be attributed to the young stellar population, but it also includes a component for the older stellar population. We have also fitted the Chary & Elbaz (2001) templates to our own de-blended SPIRE data, calculated the IR luminosities, and converted these into SFRs. From this fitting, we recover similar SFRs to those found by Schreiber et al. (2015; Fig. C.1), blue triangles. Hence, the SFRs from SED fitting are systematically lower than those found using IR luminosity conversion.

thumbnail Fig. C.1

Comparisons of Schreiber et al. (2015, red crosses) SFR-M⋆ positions with redshift to different methods of deriving SFR in this work. The mean CIGALE SFR in each mass bin is in magenta, the fitting of the de-blended SPIRE data with Chary & Elbaz (2001, CE01) templates to determine the infrared luminosity (LIR) derived SFR is indicated with blue triangles, the IAS Stacking Library derived SFR with yellow stars, the full Schreiber et al. (2015) MS trends with red dashed lines, and our own fits to the Schreiber et al. (2015) data with orange dot-dashed lines. Horizontal error bars are omitted for ease of examination, but are all 0.25 dex. The Schreiber et al. (2015) data is complete to lower masses due to their use of the deeper GOODS data. As can be seen, the SFRs from SED fitting are systematically lower than those from stacking or template fitting. However, our stacked and CE01 data points are consistent with Schreiber et al. (2015) within error.

The forward modelling was re-done using the SFRs derived from the fitting of the Chary & Elbaz (2001) templates to our SPIRE data. The slope (α) and normalisation (β) were found to evolve as (C.1)

As can be seen in Fig. C.2, the slope using the Chary & Elbaz (2001) SFR evolves with redshift at the same rate as when the CIGALE SFRs are used but the slopes are steeper, a result of high-mass galaxies having more older stars (e.g. Gallazzi et al. 2005), which will lower the CIGALE SFR with respect to the Chary & Elbaz (2001) SFR for these high-mass objects. As expected, the normalisations are closer to those of the Schreiber et al. (2015) low-mass slope and become consistent, within error, at high redshift.

The remaining normalisation difference at low redshift is a result of the enforced slope of unity in Schreiber et al. (2015). We fit Eq. (8) to the data in Schreiber et al. (2015) Fig. 10 and find slopes below unity. As a result of this fitting, the normalisations from our fits to the Schreiber et al. (2015) data agree with the normalisations we find when fitting to our IR luminosity derived SFRs. Figure C.2b also includes the normalisations for the full mass range from Schreiber et al. (2015), including their high-mass turn-over, as orange crosses. These normalisations are in good agreement with the normalisations found in this work, demonstrating how enforcing a low-mass slope of unity results in normalisations that are too high at lower redshifts.

thumbnail Fig. C.2

Comparison of the UVJ selected MS results of this work using the Chary & Elbaz (2001, CE01) template derived SFR (dark purple triangles) with the Schreiber et al. (2015) low-mass MS (red). The α and β parameters from Eq. (8) are in panels a and b, respectively, and the trends from Fig. 10are in magenta. As Schreiber et al. (2015) hold their low-mass slope constant at unity, this is indicated in panel a as the flat red line. The orange crosses in panel b are the normalisations for the full mass range (including turn-over) from Schreiber et al. (2015). The redshifts shown for this work are the mean redshift in each redshift bin while the horizontal error bars show the width of the redshift bin. The good agreement between the purple triangles and orange crosses and poor agreement between the purple triangles and red circles demonstrates how forcing a low-mass slope of unity results in normalisations that are too high at lower redshifts.

Our own fitting to the Schreiber et al. (2015) MS also provides a slightly different interpretation of their data. In addition to fitting Eq. (8) to their data, we also fit Eq. (7). The result of this refitting is that we only find evidence of a high-mass turn-over in the two lowest redshift bins used by Schreiber et al. (2015); all the other redshift bins are consistent with a simple power law (Fig. C.1 red dashed lines). Hence, their high-mass turn-over above z = 1.2 may be a result of forcing a low-mass slope of unity.

To explore whether our lack of a high-mass turn-over at low redshift is due to stacking, we stacked the SPIRE (from HerMES, Oliver et al. 2012) and PACS (from PACS Evolutionary Probe, Lutz et al. 2011) images of the COSMOS field on our source positions, using the same UVJ selection, mass bins, and redshift bins as Schreiber et al. (2015) to generate a fair comparison. We used two independent stacking codes: the IAS Stacking Library (Bavouzet 2008; Béthermin et al. 2010a) for luminosity stacking (Fig. C.1 yellow stars), and Simstack (Viero et al. 2013) for flux density stacking. The stacked Herschel flux densities were then fitted with the Chary & Elbaz (2001) IR templates to calculate the IR luminosity and SFR. We find no evidence of a high-mass turn-over using either of the stacking tools.

References

  1. Álvarez-Márquez, J., Burgarella, D., Heinis, S., et al. 2016, A&A, 587, A122 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Bavouzet, N. 2008, Theses, Université Paris Sud - Paris XI [Google Scholar]
  3. Béthermin, M., Dole, H., Beelen, A., & Aussel, H. 2010a, A&A, 512, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Béthermin, M., Dole, H., Cousin, M., & Bavouzet, N. 2010b, A&A, 516, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Boggs, P. T., & Rogers, J. E. 1990, in Statistical Analysis of Measurement Error Models and Applications (Arcata, CA, 1989) (Amer. Math. Soc., Providence, RI), Contemp. Math., 112, 183 [CrossRef] [Google Scholar]
  6. Bourne, N., Dunlop, J. S., Merlin, E., et al. 2017, MNRAS, 467, 1360 [NASA ADS] [Google Scholar]
  7. Brinchmann, J., Charlot, S., White, S. D. M., et al. 2004, MNRAS, 351, 1151 [NASA ADS] [CrossRef] [Google Scholar]
  8. Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000 [NASA ADS] [CrossRef] [Google Scholar]
  9. Burgarella, D., Buat, V., & Iglesias-Páramo, J. 2005, MNRAS, 360, 1413 [NASA ADS] [CrossRef] [Google Scholar]
  10. Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682 [NASA ADS] [CrossRef] [Google Scholar]
  11. Chabrier, G. 2003, PASP, 115, 763 [NASA ADS] [CrossRef] [Google Scholar]
  12. Charlot, S., & Fall, S. M. 2000, ApJ, 539, 718 [NASA ADS] [CrossRef] [Google Scholar]
  13. Chary, R., & Elbaz, D. 2001, ApJ, 556, 562 [NASA ADS] [CrossRef] [Google Scholar]
  14. Ciesla, L., Charmandaris, V., Georgakakis, A., et al. 2015, A&A, 576, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Dale, D. A., Cohen, S. A., Johnson, L. C., et al. 2009, ApJ, 703, 517 [NASA ADS] [CrossRef] [Google Scholar]
  16. Draine, B. T., & Li, A. 2007, ApJ, 657, 810 [NASA ADS] [CrossRef] [Google Scholar]
  17. Draine, B. T., Aniano, G., Krause, O., et al. 2014, ApJ, 780, 172 [NASA ADS] [CrossRef] [Google Scholar]
  18. Dunlop, J. S., McLure, R. J., Biggs, A. D., et al. 2017, MNRAS, 466, 861 [NASA ADS] [CrossRef] [Google Scholar]
  19. Dunne, L., Gomez, H. L., da Cunha, E., et al. 2011, MNRAS, 417, 1510 [NASA ADS] [CrossRef] [Google Scholar]
  20. Elbaz, D., Daddi, E., Le Borgne, D., et al. 2007, A&A, 468, 33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Foreman-Mackey, D. 2016, The Journal of Open Source Software, 24 [Google Scholar]
  22. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [CrossRef] [Google Scholar]
  23. Fritz, J., Franceschini, A., & Hatziminaoglou, E. 2006, MNRAS, 366, 767 [Google Scholar]
  24. Furlong, M., Bower, R. G., Theuns, T., et al. 2015, MNRAS, 450, 4486 [Google Scholar]
  25. Gallazzi, A., Charlot, S., Brinchmann, J., White, S. D. M., & Tremonti, C. A. 2005, MNRAS, 362, 41 [NASA ADS] [CrossRef] [Google Scholar]
  26. Genel, S., Vogelsberger, M., Springel, V., et al. 2014, MNRAS, 445, 175 [Google Scholar]
  27. Genzel, R., Tacconi, L. J., Lutz, D., et al. 2015, ApJ, 800, 20 [NASA ADS] [CrossRef] [Google Scholar]
  28. Górski, K. M., Hivon, E., Banday, A. J., et al. 2005, ApJ, 622, 759 [NASA ADS] [CrossRef] [Google Scholar]
  29. Griffin, M. J., Abergel, A., Abreu, A., et al. 2010, A&A, 518, L3 [Google Scholar]
  30. Harris, K., Farrah, D., Schulz, B., et al. 2016, MNRAS, 457, 4179 [Google Scholar]
  31. Hurley, P. D., Oliver, S., Betancourt, M., et al. 2017, MNRAS, 464, 885 [NASA ADS] [CrossRef] [Google Scholar]
  32. Johnston, R., Vaccari, M., Jarvis, M., et al. 2015, MNRAS, 453, 2540 [NASA ADS] [CrossRef] [Google Scholar]
  33. Jones, E., Oliphant, T., Peterson, P., et al. 2001, SciPy: Open source scientific tools for Python [Google Scholar]
  34. Kennicutt, R. C., & Evans, N. J. 2012, ARA&A, 50, 531 [NASA ADS] [CrossRef] [Google Scholar]
  35. Komatsu, E., Smith, K. M., Dunkley, J., et al. 2011, ApJS, 192, 18 [NASA ADS] [CrossRef] [Google Scholar]
  36. Kurczynski, P., & Gawiser, E. 2010, AJ, 139, 1592 [Google Scholar]
  37. Laidler, V. G., Papovich, C., Grogin, N. A., et al. 2007, PASP, 119, 1325 [Google Scholar]
  38. Laigle, C., McCracken, H. J., Ilbert, O., et al. 2016, ApJS, 224, 24 [NASA ADS] [CrossRef] [Google Scholar]
  39. Larson, D., Dunkley, J., Hinshaw, G., et al. 2011, ApJS, 192, 16 [NASA ADS] [CrossRef] [Google Scholar]
  40. Lee, N., Sanders, D. B., Casey, C. M., et al. 2015, ApJ, 801, 80 [NASA ADS] [CrossRef] [Google Scholar]
  41. Lo Faro, B., Buat, V., Roehlly, Y., et al. 2017, MNRAS, 472, 1372 [NASA ADS] [CrossRef] [Google Scholar]
  42. Lutz, D., Poglitsch, A., Altieri, B., et al. 2011, A&A, 532, A90 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Madau, P., & Dickinson, M. 2014, ARA&A, 52, 415 [NASA ADS] [CrossRef] [Google Scholar]
  44. McCracken, H. J., Milvang-Jensen, B., Dunlop, J., et al. 2012, A&A, 544, A156 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  45. Merlin, E., Fontana, A., Ferguson, H. C., et al. 2015, A&A, 582, A15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  46. Meurer, G. R., Heckman, T. M., & Calzetti, D. 1999, ApJ, 521, 64 [NASA ADS] [CrossRef] [Google Scholar]
  47. Mitchell,P. D., Lacey, C. G., Baugh, C. M., & Cole, S. 2013, MNRAS, 435, 87 [NASA ADS] [CrossRef] [Google Scholar]
  48. Mitra, S., Davé, R., Simha, V., & Finlator, K. 2017, MNRAS, 464, 2766 [NASA ADS] [CrossRef] [Google Scholar]
  49. Nguyen, H. T., Schulz, B., Levenson, L., et al. 2010, A&A, 518, L5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  50. Noeske, K. G., Weiner, B. J., Faber, S. M., et al. 2007, ApJ, 660, L43 [Google Scholar]
  51. Noll, S., Burgarella, D., Giovannoli, E., et al. 2009, A&A, 507, 1793 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Oliver, S. J., Bock, J., Altieri, B., et al. 2012, MNRAS, 424, 1614 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Pan, Z., Zheng, X., & Kong, X. 2017, ApJ, 834, 39 [NASA ADS] [CrossRef] [Google Scholar]
  54. Pannella, M., Elbaz, D., Daddi, E., et al. 2015, ApJ, 807, 141 [NASA ADS] [CrossRef] [Google Scholar]
  55. Pearson, W. J., Wang, L., van der Tak, F. F. S., et al. 2017, A&A, 603, A102 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  56. Poglitsch, A., Waelkens, C., Geis, N., et al. 2010, A&A, 518, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  57. Pozzetti, L., Bolzonella, M., Zucca, E., et al. 2010, A&A, 523, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  58. Rieke, G. H., Young, E. T., Engelbracht, C. W., et al. 2004, ApJS, 154, 25 [NASA ADS] [CrossRef] [Google Scholar]
  59. Roseboom, I. G., Oliver, S. J., Kunz, M., et al. 2010, MNRAS, 409, 48 [NASA ADS] [CrossRef] [Google Scholar]
  60. Roseboom, I. G., Ivison, R. J., Greve, T. R., et al. 2012, MNRAS, 419, 2758 [NASA ADS] [CrossRef] [Google Scholar]
  61. Safarzadeh, M., Ferguson, H. C., Lu, Y., Inami, H., & Somerville, R. S. 2015, ApJ, 798, 91 [NASA ADS] [CrossRef] [Google Scholar]
  62. Salpeter, E. E. 1955, ApJ, 121, 161 [Google Scholar]
  63. Sanders, D. B., Salvato, M., Aussel, H., et al. 2007, ApJS, 172, 86 [Google Scholar]
  64. Schreiber, C., Pannella, M., Elbaz, D., et al. 2015, A&A, 575, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Scoville, N., Aussel, H., Brusa, M., et al. 2007, ApJS, 172, 1 [NASA ADS] [CrossRef] [Google Scholar]
  66. Scoville, N., Aussel, H., Sheth, K., et al. 2014, ApJ, 783, 84 [NASA ADS] [CrossRef] [Google Scholar]
  67. Scoville, N., Sheth, K., Aussel, H., et al. 2016, ApJ, 820, 83 [NASA ADS] [CrossRef] [Google Scholar]
  68. Scudder, J. M., Oliver, S., Hurley, P. D., et al. 2016, MNRAS, 460, 1119 [NASA ADS] [CrossRef] [Google Scholar]
  69. Serra, P., Amblard, A., Temi, P., et al. 2011, ApJ, 740, 22 [NASA ADS] [CrossRef] [Google Scholar]
  70. Sparre, M., Hayward, C. C., Springel, V., et al. 2015, MNRAS, 447, 3548 [Google Scholar]
  71. Speagle, J. S., Steinhardt, C. L., Capak, P. L., & Silverman, J. D. 2014, ApJS, 214, 15 [NASA ADS] [CrossRef] [Google Scholar]
  72. Tacchella, S., Dekel, A., Carollo, C. M., et al. 2016, MNRAS, 457, 2790 [NASA ADS] [CrossRef] [Google Scholar]
  73. Tacconi, L. J., Genzel, R., Neri, R., et al. 2010, Nature, 463, 781 [Google Scholar]
  74. Tomczak, A. R., Quadri, R. F., Tran, K.-V. H., et al. 2016, ApJ, 817, 118 [NASA ADS] [CrossRef] [Google Scholar]
  75. Viero, M. P., Moncelsi, L., Quadri, R. F., et al. 2013, ApJ, 779, 32 [NASA ADS] [CrossRef] [Google Scholar]
  76. Vogelsberger, M., Genel, S., Springel, V., et al. 2014, MNRAS, 444, 1518 [NASA ADS] [CrossRef] [Google Scholar]
  77. Wang, L., Viero, M., Clarke, C., et al. 2014, MNRAS, 444, 2870 [NASA ADS] [CrossRef] [Google Scholar]
  78. Wang, L., Norberg, P., Bethermin, M., et al. 2016, A&A, 592, L5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  79. Whitaker, K. E., Labbé, I., van Dokkum, P. G., et al. 2011, ApJ, 735, 86 [NASA ADS] [CrossRef] [Google Scholar]
  80. Whitaker, K. E., van Dokkum, P. G., Brammer, G., & Franx, M. 2012, ApJ, 754, L29 [NASA ADS] [CrossRef] [Google Scholar]
  81. Wright, A. H., Robotham, A. S. G., Bourne, N., et al. 2016, MNRAS, 460, 765 [NASA ADS] [CrossRef] [Google Scholar]
1

It is possible to recover the scatter of the objects in a stack. However, this requires strong assumptions on the distribution of the galaxies in the stack (Schreiber et al. 2015; Wang et al. 2016).

All Tables

Table 1

Telescopes and associated bands that were used for CIGALE spectral energy distribution fitting from the COSMOS2015 catalogue.

In the text
Table 2

Parameters from fitting Eq. (8) to the UVJ selected star forming galaxies, with α the slope and β the normalisation at log(M⋆/M⊙) = 10.5.

In the text
Table A.1

Parameters used for the various properties in the CIGALE model SEDs where they differ from the default values.

In the text

All Figures

thumbnail Fig. 1

Imageof the SPIRE 250 μm COSMOS coverage (Oliver et al. 2012, 7.84 deg2, red) with the overlayed MIPS 24 μm COSMOS coverage (Sanders et al. 2007, 2.23 deg2, blue). The data used in this work were cut to match the MIPS 24 μm coverage.
In the text
thumbnail Fig. 2

Brief summary of the science pipeline.
In the text
thumbnail Fig. 3

Scatter of the 250 μm residual map using different depth cuts on the prior list. The blue and red points are with and without 3σ clipping, respectively, while the green line is the 1σ instrument noise for the COSMOS field.
In the text
thumbnail Fig. 4

Masses of all objects detected in the Ks band (blue) are shown against redshift along with the faintest 20% in each redshift bin (red) for the (panel a) deep and (panel b) ultra-deep regions. The 90% completeness limit is shown by the thick black lines, while the dashed black lines show the edges of each redshift bin.
In the text
thumbnail Fig. 5

Example of the data generated at one step of the MCMC routine for the lowest (0.2 ≤ z < 0.5) redshift bin, shown as number density from high (dark red) to low (dark blue). The MS being tested at this step, in this case the most likely step, is shown as the red line, while the contours for number density of the observed data are shown in orange. The size of the average observed error on SFR and M⋆ is also shown as a blue cross.
In the text
thumbnail Fig. 6

Corner plot (Foreman-Mackey 2016) of the marginalised posterior of the forward modellingroutine applied to a mock data set with known slope (0.6), normalisation (0.7 log (M⊙∕yr−1)), and scatter (0.3 dex), shown as the blue lines. Panels on the diagonal show the 1D marginalised posteriors for the slope, normalisation and scatter (left to right panels). Off-diagonal panels show the combined 2D posteriors as labelled by their axes. The recovered 16th, 50th, and 84th percentiles are shown by the dashed vertical lines; all the input parameters are recovered, within error.
In the text
thumbnail Fig. 7

Comparison of fitting with Eq. (7) (orange dashed line) and Eq. (8) (red line) in the 0.5 ≤ z < 0.8 redshift bin. The galaxies are shown as a number density plot, with dark red being high number density and dark blue low number density, and the size of the average error on SFR and M⋆ is shown as a blue cross. There is very little difference in shape between the two fits, demonstrating that Eq. (8) is the preferred form of the MS.
In the text
thumbnail Fig. 8

Fitting of Eq. (8) to the objects in the redshift bins as labelled. The solid line is the most likely MS across the fitted M⋆ range, while the dotted line is an extrapolation across the M⋆ range of theplot. The galaxies are shown as a number density plot, with dark red being high number density and dark blue low number density. The vertical density discontinuities are the result of the two depths of data used: the deep mass limit is the yellow dashed line and the ultra-deep mass limit is the dot-dashed yellow line. Each panel also indicates the χ2 of the most likely MS () and the number of degrees of freedom in the fitting (ndof), and shows the size of the average SFR and M⋆ errors as a blue cross.
In the text
thumbnail Fig. 9

Fitted MS trends using Eq. (8). The solid and dashed lines are the MS across the fittedM⋆ range. The dotted lines are extrapolations across the M⋆ range of the plot. The normalisation of the MS clearly increases as redshift increases. The slight decrease in slope at low redshift can be seen along with the increase above z = 1.1. The very steep slopes at 1.8 ≤ z < 2.9 can also be seen.
In the text
thumbnail Fig. 10

Comparison of the UVJ selected MS results of this work (magenta) with the observational MS from Speagle et al. (2014, FUV data in blue; IR data in green; and radio, hydrogen lines, and UV SED fitting in yellow), Schreiber et al. (2015) low mass MS (red), and the MS from the Illustris Simulation (Vogelsberger et al. 2014; Sparre et al. 2015, black diamonds). The α and β parameters from Eq. (8) are in panels a and b, respectively. As Schreiber et al. (2015) hold their low-mass slope constant at unity, this is indicated in panel a as a flat red line. The redshifts shown for this work are the mean redshift in each redshift bin, while the horizontal error bars show the width of the redshift bin. A version of this plot using SFRs derived from IR template fitting can be found in Appendix C.
In the text
thumbnail Fig. 11

Intrinsic scatter about the MS found from the MCMC fitting of the data in each redshift bin to Eq. (8). This work (magenta) is shown along with the observed scatters from Speagle et al. (2014, FUV data in blue, IR data in green, and other data in yellow) and Schreiber et al. (2015, red) as well as the median scatter found at each redshift in the Illustris Simulation (Genel et al. 2014; Vogelsberger et al. 2014; Sparre et al. 2015, black diamonds). The redshifts shown for this work are the mean redshift in each redshift bin while the horizontal error bars show the width of the redshift bin. Also shown is the best fit to this work’s intrinsic scatter assuming no redshift evolution (dashed magenta line). The intrinsic scatter found in this work is consistent with existing literature, although above z ≈ 1.8 our intrinsicscatter is smaller than is found in works that use IR SFR traces.
In the text
thumbnail Fig. B.1

Distributions for the difference between the observed ALMA 870 μm flux densities and the predicted 870 μm flux densities from CIGALE using an error expansion factor of 2 (panel a), 3 (panel b), and 4 (panel c) in XID+. The Gaussian distribution for each expansion factor is also shown in orange, along with the means (μ) and standard deviations (σ) of the distributions. All the distributions have an approximately consistent σ so the best expansion factor was deemed to be that with μ closest to zero: 2 times the error from CIGALE.
In the text
thumbnail Fig. C.1

Comparisons of Schreiber et al. (2015, red crosses) SFR-M⋆ positions with redshift to different methods of deriving SFR in this work. The mean CIGALE SFR in each mass bin is in magenta, the fitting of the de-blended SPIRE data with Chary & Elbaz (2001, CE01) templates to determine the infrared luminosity (LIR) derived SFR is indicated with blue triangles, the IAS Stacking Library derived SFR with yellow stars, the full Schreiber et al. (2015) MS trends with red dashed lines, and our own fits to the Schreiber et al. (2015) data with orange dot-dashed lines. Horizontal error bars are omitted for ease of examination, but are all 0.25 dex. The Schreiber et al. (2015) data is complete to lower masses due to their use of the deeper GOODS data. As can be seen, the SFRs from SED fitting are systematically lower than those from stacking or template fitting. However, our stacked and CE01 data points are consistent with Schreiber et al. (2015) within error.
In the text
thumbnail Fig. C.2

Comparison of the UVJ selected MS results of this work using the Chary & Elbaz (2001, CE01) template derived SFR (dark purple triangles) with the Schreiber et al. (2015) low-mass MS (red). The α and β parameters from Eq. (8) are in panels a and b, respectively, and the trends from Fig. 10are in magenta. As Schreiber et al. (2015) hold their low-mass slope constant at unity, this is indicated in panel a as the flat red line. The orange crosses in panel b are the normalisations for the full mass range (including turn-over) from Schreiber et al. (2015). The redshifts shown for this work are the mean redshift in each redshift bin while the horizontal error bars show the width of the redshift bin. The good agreement between the purple triangles and orange crosses and poor agreement between the purple triangles and red circles demonstrates how forcing a low-mass slope of unity results in normalisations that are too high at lower redshifts.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.