© ESO, 2017

1. Introduction

The morphology of galaxies is a visible outcome of all processes at work when galaxies are assembling along cosmic time (e.g., Conselice 2014). While at lower redshifts galaxies are predominantly regular and symmetric in the form of elliptical or spiral galaxies, the situation at redshifts higher than the peak in star formation z ≃ 1.5–2 is very different with galaxies mostly having irregular shapes (e.g., Buitrago et al. 2013; Mortlock et al. 2013; Huertas-Company et al. 2015).

Quantitative morphology allows us to follow the evolution of the two main components of galaxies, disks and bulges, and studies have focused on understanding the evolutionary paths of these fundamental components (e.g., Bruce et al. 2012, 2014a,b; Lang et al. 2014; Tasca et al. 2014b; Margalef-Bentabol et al. 2016). Observational studies and numerical simulations seem to support a picture where the original baryonic matter collapsing in a dark matter halo conserves its angular momentum which leads to the formation a rotation supported proto-disk. In this scenario, violent disk instabilities (VDI) develop and then form large star-forming clumps. These large clumps may then migrate to the disk center to form the proto-bulges of spiral galaxies (e.g., Genzel et al. 2008; Elmegreen et al. 2008, 2009; Dekel & Krumholz 2013; Perez et al. 2013; Bournaud et al. 2014; Bournaud 2016). Alternatively, considerable effort has been invested in studying the development of elliptical galaxies along cosmic time. On this scenario, bulges would form from violent dissipative events, going through a stage of compact quiescent objects at z ~ 2. Such objects would then grow into massive ellipticals after star formation quenching and a succession of merging events (e.g., Hopkins et al. 2008; Oser et al. 2012; López-Sanjuan et al. 2012; Porter et al. 2014). From a theory and simulations perspective, two major physical processes are therefore thought to play a key role in the evolution of different types of galaxies: gas accretion mostly in cold form and flowing along the cosmic web and/or smoothly distributed (e.g., Dekel et al. 2009; Rodriguez-Gomez et al. 2016), and major or minor merging (e.g., Hopkins et al. 2008; Di Matteo et al. 2008), each leading to different morphological signatures.

An important morphological property related to these processes is the number of large clumps of gas and stars present within a galaxy. Clumps may have either an in-situ or an ex situ formation history depending on which physical processes are at play (e.g., Shibuya et al. 2016). Clumps forming from disk instabilities are mostly the result of in-situ dynamical processes, while major and minor merging events are driving clumps of ex situ material into a galaxy in assembly phase. The properties and fate of these clumps have been extensively studied in the literature both from observations (e.g., Elmegreen et al. 2004; Elmegreen & Elmegreen 2006; Elmegreen et al. 2007, 2013; Genzel et al. 2008; Bournaud et al. 2008; Wuyts et al. 2012; Guo et al. 2012) and in hydrodynamical simulations (e.g., Bournaud et al. 2007; Dekel et al. 2009; Genel et al. 2012; Bournaud et al. 2014; Mandelker et al. 2014; Tamburello et al. 2015; Oklopčić et al. 2017). Both report clumps with typical stellar masses of 10⁷−10⁹ M_⊙ and with short to long lifetimes of ≲1 Gyr. Some studies suggest that short-lived clumps are destroyed by internal feedback before migrating to the center of the galaxy (e.g., Genel et al. 2012).

The measurement of clumps at redshifts beyond the peak of cosmic star formation is notoriously difficult. The identification of clumpy galaxies in the dominant population of irregular galaxies at high redshift dates back to the first deep Hubble Space Telescope (HST) images (e.g., Williams et al. 1996). Subsequent studies at z ~ 1−3 revealed that clumpy galaxies are more numerous than in the local Universe (e.g., Abraham et al. 1996; van den Bergh et al. 1996; Giavalisco et al. 1996; Elmegreen et al. 2004, 2007, 2008, 2009; Elmegreen & Elmegreen 2006; Elmegreen et al. 2013; Kubo et al. 2013, 2016; Glazebrook 2013; Tadaki et al. 2014; Murata et al. 2014; Guo et al. 2015; Garland et al. 2015; Bournaud 2016). The fraction of clumpy galaxies was studied from z = 0 up to the most distant galaxies identified today at z ~ 10 (Ravindranath et al. 2006; Guo et al. 2012, 2015; Shibuya et al. 2016). In the sample of galaxies with 0.5 <z < 3 in the CANDELS survey, Guo et al. (2015) find that galaxies with low mass log (M_⋆/M_⊙) < 9.8 have a high fraction of off-center clumps f_clumpy ~ 60% constant over the observed redshift range, while for intermediate and massive galaxies f_clumpy decreases from ~40% at z ~ 3 to ~15% at z ~ 0.5. Combining deep HST imaging from the 3D-HST, CANDELS, HUDF and HFF surveys Shibuya et al. (2016) claim that ${\mathit{f}}_{\mathrm{clumpy}}^{\mathrm{UV}}$ follows an evolution similar to that of the star formation rate density (e.g., Madau & Dickinson 2014), increasing from z ≃ 8 to z ≃ 1−3 and decreasing from z ≃ 1 to z ≃ 0. In comparing these observed trends with the predictions of simulations, Guo et al. (2015) conclude that VDI are likely responsible for f_clumpy at high mass and that minor mergers are a viable explanation for f_clumpy at intermediate mass for z < 1.5, while Shibuya et al. (2016) conclude that VDI is the main origin of all clumps. Both these studies exclude major mergers as a possible contribution to the evolution of the fraction of clumps. This is somewhat surprising as the major merger fraction for star-forming galaxies is observed to increase up to ~20% at z ≃ 1−2 (e.g., Lotz et al. 2011; López-Sanjuan et al. 2011, 2013), staying high at least to z ~ 3−4 (e.g., Conselice et al. 2008; Tasca et al. 2014a), and one would therefore expect that a fraction of the clumps is related to ex situ major merging.

When deriving volume average quantities like f_clumpy it is critical to observe samples representative of the general star-forming population, with well defined physical parameters including luminosity, stellar masses and sizes, and following a well defined selection function. The size of irregular high redshift galaxies is very different when using their effective radius or a total size (Ribeiro et al. 2016), and searching for clumps using one or the other of these sizes to define a search area may lead to significantly different results on f_clumpy. We show that total sizes are larger by a factor of ~1.5−2 than what previous findings using r_e claim, and that this difference is larger at higher redshifts (comparing z ~ 2 to z ~ 4, Ribeiro et al. 2016). Moreover, surface brightness dimming and luminosity evolution with redshift complicate the issue of the definition and identification of a clump. One must be aware of the impact of specific definitions and sample selection, for subsequent analysis to have a statistically robust meaning in a representative volume of the Universe at a given redshift.

In this paper we have used a sample of 1242 galaxies with secure spectroscopic redshifts from the VIMOS Ultra-Deep Survey (VUDS) and HST imaging in the COSMOS and ECDFS fields to measure f_clumpy. We used a clump definition expanded from Guo et al. (2015), counting all clumps within the isophote defining the total size r_T,100 and in a 20 kpc physical radius. In our clump definition we made no attempt at distinguishing between galaxy nucleus and a star-forming clump. The reasoning behind this choice of classification is that at these redshifts the nucleus of one galaxy may well be in itself a star-forming clump, as the timescales involved in clump formation and respective lifetimes can be comparable to the age of the galaxy we are observing. We then used the observed properties of the detected star-forming regions to infer their possible nature and formation mechanisms.

This paper is organized as follows. The sample and the data we used are described in Sect. 2. We present our method for identifying and characterizing clumps in Sect. 3. The evolution of the fraction of clumps is presented in Sect. 4. We analyze the statistical properties of clumps in star-forming galaxies in Sect. 5 where we explore the distribution of clumps in galaxies, their areas and luminosities and their spatial distribution. The results are discussed in Sect. 6 and summarized in Sect. 7. We use a cosmology with H₀ = 70 km s^-1 Mpc^-1, Ω_0,Λ = 0.7 and Ω_0,m = 0.3. All magnitudes are given in the AB system.

2. Data and sample selection

We aimed to work with a sample of galaxies with accurate spectroscopic redshifts. This avoids the uncertainties generally related to photometric redshifts with typical accuracy at z ~ 3.5 of δz ~ 0.2−0.3 when deriving important physical quantities like sizes, the luminosity of an isophote, or total luminosities and stellar masses.

We then imposed that the galaxies with spectroscopic redshifts have deep HST imaging available. As algorithms for clump detection are sensitive to the redshift of the galaxy (see Sect. 3), the HST images combined with the spectroscopic redshift of a galaxy enable the computation of reliable isophotal areas within which we identified clumps and quantify their properties. We defined as a clump any component in a galaxy image that counts more than five connected pixels above a reference isophote (see Sect. 3).

We further restricted our sample above a given stellar mass defined as a function of redshift from the evolving mass function. This ensures that we followed a similar population of galaxies over the explored redshift range. This selection of the sample and associated data is described in detail in Sects. 2.1 to 2.3.

2.1. Starting sample selection: galaxies with spectroscopic redshifts from the VIMOS Ultra Deep Survey (VUDS)

We use the VUDS sample, the spectroscopic survey conducted by our team with ~10 000 objects targeted in an area of 1 deg² on three separate fields: COSMOS, ECDFS and VVDS-02h (Le Fèvre et al. 2015). The VUDS sample covers the redshift range 2 <z < 6 +, and produces a magnitude-selected sample of galaxies with 22.5 ≤ i_AB ≤ 25 and spectroscopic redshifts.

The spectra were obtained using the VIMOS spectrograph on the ESO-VLT (Le Fèvre et al. 2003) covering, with two low resolution grisms (R = 230), a wavelength range of 3650 Å <λ < 9350 Å. The total integration time is T_exp ~ 14 h per pointing and grism. Data processing was performed within the VIPGI environment (Scodeggio et al. 2005), and followed by extensive spectroscopic redshift measurements campaigns using the EZ redshift measurement engine (Garilli et al. 2010). At the end of this process each galaxy has a spectroscopic redshift measurement, and an associated redshift reliability flag.

In the following we analyze exclusively VUDS galaxies with the most reliable spectroscopic redshifts, that is, galaxies with spectroscopic redshift flags 2 and 9 (~75% probability to be correct), and flags 3 and 4 (95 to 100% probability to be correct). For more information on this sample we refer the reader to Le Fèvre et al. (2015).

2.2. Sample with imaging data in the VUDS-COSMOS and VUDS-ECDFS fields

To reliably identify clumps in galaxies we use those galaxies in VUDS as defined in Sect. 2.1 which have deep high resolution images available from HST surveys. The COSMOS survey (Scoville et al. 2007; Koekemoer et al. 2007) has ≈2deg² of HST ACS F814W (I band) coverage down to a limiting magnitude of 27.2 (5σ point-source detection limit), and covers the entirety of the VUDS COSMOS area of ~1800 square arcmin (Le Fèvre et al. 2015). This field provides the largest HST imaging data set of galaxies with spectroscopic redshifts. In the ECDFS we use the F814W images overlapping with the VUDS footprint over 675 square arcmin (see Le Fèvre et al. 2015). The F814W images cover the UV rest-frame from ~2700 Å at z ~ 2 to ~1500 Å at z ~ 4.5.

To probe the rest-frame optical morphology we use the smaller sample of galaxies in VUDS covered by the CANDELS HST imaging survey as presented in Tasca et al. (2017), a subset of the full CANDELS area (Grogin et al. 2011; Koekemoer et al. 2011). Even if small, the area of ~270 square arcminutes is important because the near-infrared coverage (with the F125W and F160W bands) allows us to compare optical rest-frame properties to UV rest-frame derived from the F814W imaging. Data in CANDELS reaches depths of 5σ of 28.4, 27.0 and 26.9 at F814W, F125W and F160W, respectively. The typical spatial resolution of these images ranges from 0.09′′ (0.8,0.6 kpc, at z = 2.0,4.5) for the F814W band up to 0.18′′ (1.5,1.2 kpc, at z = 2.0,4.5) in the F160W band.

Fig. 1 Stellar masses and star formation rates for the galaxies in the selected sample. In each panel, the light gray points show all VUDS galaxies with measured redshift in each bin. The colored circles show the 1242 galaxies in the selected stellar mass range and with correct redshift probability of >75% (see text for more details).

2.3. Selecting a stellar-mass limited sub-sample

In the final step of our sample selection we define a mass-selected sample from the sample defined in Sect. 2.2, based on stellar masses computed from spectral energy distribution (SED) fitting on the large set of photometric data available in the COSMOS and ECDFS fields. The SED fitting procedure closely follows the method described in Ilbert et al. (2013) using the latest version of the Le Phare code (maintained by Olivier Ilbert), with the specifics for VUDS detailed in Tasca et al. (2015) and summarized in Table 1. The main parameter of interest in this paper is the stellar mass, M_⋆ for which the median values of the probability distribution function are used. We refer the reader to Ilbert et al. (2013) and Thomas et al. (2017) for typical uncertainties on these quantities (see also Tasca et al. 2015). The statistical uncertainty on stellar masses is typically ~0.1 dex. Systematic uncertainties due, for example, to the choice of initial mass function (IMF) can be as high as a factor of ~1.7, if a Salpeter (1955) IMF is used instead of the Chabrier (2003) IMF. The stellar mass selection (discussed below), color correction detailed in Sect. 3.2 and the k-correction used in Sect. 3.3 rely on the parameters derived from the SED fitting with Le Phare. Using stellar masses as derived from another SED fitting code as for example, the GOSSIP+ code (see Thomas et al. 2017), performing SED fitting on the joint photometry and spectroscopy data, has no impact on the results presented in this paper.

We present the stellar mass–star forming rate (SFR) relation of our sample divided in four redshift bins in Fig. 1. It shows that our sample is representative of the high stellar mass population in VUDS and that it probes galaxies with typical median stellar masses of 10¹⁰ M_⊙ and ranging from log ₁₀(M_⋆/M_⊙) ≳ 9.35 and up to log ₁₀(M_⋆/M_⊙) ≳ 11. In terms of SFR, our galaxies are in the range 0.5 ≲ log ₁₀(SFR) ≲ 3.0 at all redshifts with median values of log ₁₀(SFR) ~ 1.4.

To follow the evolution of galaxy properties in a similar population as a function of redshift, broadly representing the same coeval population, we defined our sample imposing an evolving lower stellar mass limit set to follow the general stellar mass growth of star-forming galaxies. We defined the lower stellar mass limit in our sample as the limit below which the VUDS sample is becoming incomplete at z = 4.5, setting log ₁₀(M_⋆/M_⊙) > 9.35. We then used the stellar mass function evolution from Ilbert et al. (2013) together with the typical sSFR of VUDS galaxies (Tasca et al. 2015) to follow the typical stellar mass growth of VUDS galaxies, and defined the stellar mass selection threshold at different redshifts using ${\log }_{\mathrm{10}}\mathrm{\left(}{\mathit{M}}_{\mathit{\star }}\mathit{/}{\mathit{M}}_{\mathrm{\odot }}\mathrm{\right)}\mathit{>}\mathrm{-}\mathrm{0.204}\mathrm{\times }\mathrm{(}\mathit{z}\mathrm{-}\mathrm{4.5}\mathrm{)}\mathrm{+}\mathrm{9.35.}$(1)We refer to Ribeiro et al. (2016) for more details on the choice of this mass selection in the context of the VUDS survey. In the higher redshift interval z > 4.5 we selected all galaxies with spectroscopic redshift and stellar masses log ₁₀(M_⋆/M_⊙) > 9 to build a sufficiently large sample to obtain a measurement of the clumpy fraction and clump properties, even if this sample is slightly incomplete.

2.4. Final sample

We identified a total of 1242 galaxies with 2 < z < 4.5 that satisfy the mass-selection and the redshift quality criteria in the VUDS survey with ACS F814W images, of which 1087 are in the COSMOS field and 155 are in the CANDELS-ECDFS field. We also examined a higher redshift sample of 96 galaxies with 4.5 < z < 6, even if the available imaging data for this higher redshift sample is not deep enough to search for clumpy galaxies at the same luminosity limit than the lower redshift sample, and we take this into consideration in our analysis.

Fig. 2 Exemplification of the method for detecting galaxy clumps. On the left we show the number of disconnected regions with more than ten pixels as a function of flux level cut. When f = 1 the isophote limit represents the limiting area of the galaxy. On the right we show the contours corresponding to three different flux level cuts at 20%, 50% and 100% of total flux that are identified in the same color as the vertical lines in the left panel.

3. Defining clumps in galaxies

There are various methods in the literature that have the goal of detecting and characterizing galactic clumps at high redshift ranging from visual classification to automated algorithms (e.g., Elmegreen et al. 2007; Puech 2010; Guo et al. 2012, 2015; Wuyts et al. 2012; Tadaki et al. 2014; Shibuya et al. 2016). The clump detection algorithms generally involve the subtraction of a smoothed version of an image from the original and then search for positive detections in the residual image that can be considered as galactic clumps, a simple and robust method. In the work by Guo et al. (2015) and Shibuya et al. (2016) a clump is detected if five connected pixels are detected at 2σ above the observed image background. A surface brightness limit defined with respect to the background corresponds to a fainter intrinsic (absolute) luminosity at lower redshift than at higher redshift because of surface brightness dimming varying as (1 + z)^-3 (Ribeiro et al. 2016). For instance, between z = 2 and z = 4 using the same surface brightness limit defined with respect to the background leads to almost a factor of five difference in the intrinsic surface brightness of these galaxies. Counting clumps may therefore be significantly affected depending on the definition of the limiting surface brightness limit. The overall lowering of the clump fraction with increasing redshift as reported by previous studies may simply be due to the bias produced by using the same observed surface brightness at all redshifts when identifying clumps. Because of the limits in detection at the highest redshifts one probes only the brighter and rarer clumps, while at the same limiting surface brightness one finds fainter clumps which are more numerous, hence mimicking a clump fraction decrease with increasing redshift.

To follow on from these considerations, we followed a slightly different approach in this paper than in the recent studies of Guo et al. (2015) and Shibuya et al. (2016). We imposed a surface brightness limit when detecting clumps which has the same intrinsic luminosity for all measured galaxies, whatever their redshift. We have therefore taken into account surface brightness dimming differences, which are important over the redshift range considered. Furthermore, the average luminosity of a galaxy population changes with redshift as shown by luminosity function studies (and the evolution of the luminosity density), and we further adjusted the limiting surface brightness when searching for clumps taking into account this redshift dependent luminosity evolution.

In this section we describe the method we follow to identify clumps and measure their luminosity and size. We also reproduced the method used by Guo et al. (2015) and Shibuya et al. (2016) for comparison with our results on the same galaxy sample.

3.1. A clump finding method consistent at all redshifts

We summarize our method in the logical steps below.

• 1.

  We smoothed a 6″ × 6″ image centered on the VUDS galaxy with a Gaussian kernel of σ = 1 pixel (0.03′′).

• 2.

  We computed the surface brightness threshold corresponding to the galaxy, from a reference intrinsic threshold imposed on all galaxies, that takes into account the relative surface brightness diming as well as the average luminosity evolution of the galaxy population.

• 3.

  We created a segmentation map that encodes which pixels are found to be above the defined surface brightness threshold.

• 4.

  From the segmentation map we selected the area of connected pixels above the adopted surface brightness threshold which contains the brightest pixel within a 0.5 arcsec radius from the target coordinates.

• 5.

  Using the selected pixels, we computed the total flux of the galaxy by summing up the flux in all the pixels.

• 6.

  We iterated (using 1000 steps) a flux level f_i from the maximum to the minimum flux detected. At each step we computed the number of disconnected regions (defined from all pixels with fluxes greater than f_i) within the segmentation map and stored the position of the local maximum associated with each disconnected region (see Fig. 2).

• 7.

  We defined the number of clumps as the number of disconnected regions detected above the isophote level corresponding to 80% of the total galaxy flux, which are associated with the different local maxima.

We repeated this process for all galaxies found in the image stamp. We note that we have only considered clumps with a detection area greater that ten pixels. This minimum area for clump detection ensures that we have excluded the contribution from random fluctuations of photon counts and also imposes a lower limit of ~3 pixels (roughly the FWHM of the PSF) on the separation of any given two clumps that we can detect.

We describe some of the key elements of this algorithm in Sects. 3.2 to 3.6.

3.2. Limiting surface brightness: isophotal threshold

One of the key elements of this algorithm is the limiting isophote at which we define the total extent of a single galaxy. For our method, it is defined as $\mathit{k}\mathrm{=}{\mathit{k}}_{\mathrm{p}}{\left(\frac{\mathrm{1}\mathrm{+}\mathit{z}}{\mathrm{1}\mathrm{+}{\mathit{z}}_{\mathrm{p}}}\right)}^{-3}\mathrm{\times }{\mathrm{10}}^{\mathrm{-}\mathrm{0.4}\mathrm{(}{\mathit{I}}_{\mathrm{obs}}\mathrm{-}{\mathit{B}}_{\mathrm{rest}}\mathrm{)}}\mathrm{\times }\frac{\mathit{L}\mathrm{\left(}\mathit{z}\mathrm{\right)}}{\mathit{L}\mathrm{\left(}{\mathit{z}}_{\mathrm{p}}\mathrm{\right)}}\mathit{,}$(2)where k_p (pivot threshold) is the value of k at z = z_p (pivot redshift, with z_p = 2 in this work) and I_obs−B_rest = 0, z is the spectroscopic redshift, I_obs is the apparent magnitude of the galaxy, B_rest is the apparent rest frame magnitude of the band B, and L(z) is the typical luminosity at redshift z (see Eq. (3) for more details). Equation (2) is composed of three different terms:

The dimming correction term, [(1 + z)/(1 + z_p)] ^-3, which accounts for the fact that due to the expansion of the Universe the observed number of photons of a source with constant intrinsic luminosity scales with (1 + z)^-3 (e.g., Ribeiro et al. 2016).

The color correction term, 10^{− 0.4(I_obs−B_rest)}, which accounts for the fact that we are not observing galaxies exactly at the same rest-frame wavelengths and corrects for the relative brightness based on the slope of the SED close to the band we are observing. To be consistent in our large redshift range we chose to use rest-frame far ultraviolet (FUV) magnitudes to compute the color correction term since it is the central wavelength at the peak of the redshift distribution of our sample.

The luminosity correction term, L(z) /L(z_p), which accounts for the average luminosity evolution in galaxies. L(z) is derived from the L^∗ values computed in FUV by Reddy & Steidel (2009) at z < 3 and by Bouwens et al. (2015) at z > 4. We defined a luminosity evolution which is steeper below z < 3 and flattens at z > 3: $\mathit{L}\mathrm{\left(}\mathit{z}\mathrm{\right)}\mathrm{=}\{\begin{array}{c}\\ & \mathit{z}\mathrm{\le }\mathrm{3}\\ & \mathit{z}\mathit{>}\mathrm{3.}\end{array}$(3)In Eq. (2) there are two parameters that are correlated. The value k_p represents the detection threshold, in units of σ, at z = z_p. In this study we set z_p = 2 and used several values of k_p throughout the paper. In particular we used k_p = 1.5,2.0 for 2 <z < 4.5 and k_p = 5.0 for 4.5 <z < 6. Since the value of k can vary by a factor of ((1 + 2)/(1 + 6))^-3 ~ 0.08 from z = 2 to z = 6 based only on surface brightness dimming we chose different values of k_p to allow for an estimate of the total extent of galaxies at the lower limit while preventing noise contaminated detections at the high redshift end.

3.3. Luminosity and sizes

In order to obtain further information on the physical characteristics of clumps we adapted the method defined by Freeman et al. (2013) to compute the intensity parameter I giving the luminosity of each clump. This algorithm defines groups of pixels associated with a local maximum inside a segmentation map. First, the image was smoothed via a Gaussian kernel with a given σ (we use σ = 1 pixel). Then, for each pixel in the segmentation map, we computed the corresponding local maximum via maximum gradient paths. This means that for each pixel we computed the intensities of the eight surrounding pixels and then moved to the pixel with maximum intensity. We continued this process until a local maximum was reached and none of the eight neighbor pixels had higher intensities than the one we were considering. So, a pixel group is defined as the collection of pixels linked to a particular local maximum.

Once we had all the pixel groups defined, we matched the local maximum of each group to the local maximum of the detected clumps as described in Sect. 3. Then we computed the luminosity and area associated with each clump by counting all the light and pixels in each group, respectively. It is important to stress that since we stopped at the border of two contiguous groups of pixels, the derived areas and luminosities are lower limits to the real values. This happens since we only considered flux from pixels where the contribution from the clump being considered is dominant with respect to neighboring clumps.

With the knowledge of the spectroscopic redshift from VUDS, we were then able compute physical areas in kpc². Absolute magnitudes are computed using the flux measured from the F814W images. We corrected the magnitudes to FUV rest-frame using the k-correction derived for the galaxy from SED fitting. This means that the clump magnitude is derived as ${\mathit{M}}_{\mathrm{clump}\mathit{,}\mathrm{FUV}}\mathrm{=}{\mathit{M}}_{\mathrm{clump}\mathit{,}\mathit{F}\mathit{814}\mathit{W}}\mathrm{+}{\mathit{k}}_{\mathrm{corr}\mathit{,}\mathrm{FUV}}\mathit{.}$(4)

3.4. Resulting clumps and pair identification

At the end of this process a list of galaxies was identified within a 20 kpc radius, a value commonly used as two equal masses bound by gravity would likely merge in a dynamical timescale of ~1 Gyr (López-Sanjuan et al. 2013; Tasca et al. 2014a), and a list of clumps was identified in each galaxy. Our algorithm may identify several galaxies within the search radius of 20 kpc, as pixels groupings isolated from other pixel groupings, that is, these groups of pixels are not connected to each other at the isophote considered. It is well possible that if deeper images were available, allowing us to lower the search isophote, some of these pixel groupings could be connected and possibly associated to a single galaxy.

To avoid ambiguous interpretation, all our results on galactic clumps take into account only clumps that are associated with the galaxy that was the target of the VUDS spectroscopy. We have spectroscopic information (important to compute k in Eq. (2)) only for these sources as well as physical parameters from SED fitting only for these galaxies. However we keep information about other galaxies detected within 20 kpc to compute a pair fraction at the same redshifts, as discussed in Sect. 6. Several examples of clumps detected in VUDS galaxies are shown in Figs. 3−5.

To transform our pair fraction into a merger fraction we corrected the observed pair count for line of sight contamination using a simple statistical argument. Starting with the expected number of galaxies within a square degree, as derived by Capak et al. (2007), we computed the probability of a galaxy of magnitude m to be found at a distance r of the VUDS target. The total number of galaxies brighter than magnitude m in one square degree was computed by integrating the number counts, N( <m) computed by Capak et al. (2007). We then estimated how many galaxies we expected to find within π × r². This is summarized by $\mathit{p}\mathrm{=}\frac{\mathit{\pi }{\mathit{r}}^{\mathrm{2}}}{\mathrm{3600}\mathrm{\times }\mathrm{3600}}\mathrm{\times }{\mathrm{\int }}_{\mathrm{15}}^{\mathit{m}}\mathit{N}\mathrm{(}\mathit{<}{\mathit{m}}^{\mathrm{\prime }}\mathrm{)}\mathrm{d}{\mathit{m}}^{\mathrm{\prime }}\mathit{.}$(5)To correct for line of sight contamination, each galaxy of magnitude m is weighted by its probability of being a contaminant, 0 <p < 1, when we computed the pair fraction.

Fig. 3 Twelve examples of galaxies for which only one clump was detected. The white circle shows the 20 kpc projected search radius. The contours represent the detection isophote of each galaxy identified within this search radius. Each circle shows the position of a detected clump within the isophote of the same color. VUDS targets are identified with red contours and circles.

Fig. 4 Twelve examples of galaxies with two clumps detected. The top six examples show cases for which at least one bright clump was detected whereas the bottom six examples show pairs of faint clumps. The white circle shows the 20 kpc projected search radius. The contours represent the detection isophote of each galaxy identified within this search radius. Each circle shows the position of a detected clump within the isophote of the same color. VUDS targets are identified with red contours and circles.

Fig. 5 Twelve examples of galaxies with more than two clumps detected. The white circle shows the 20 kpc projected search radius. The contours represent the detection isophote of each galaxy identified within this search radius. Each circle shows the position of a detected clump within the isophote of the same color. VUDS targets are identified with red contours and circles.

3.5. False clump detections

To estimate the impact of false detections on the inferred distributions of areas and luminosities we performed a simple test. We ran our clump detection algorithm using the same threshold as before on the negative image of each galaxy and computed the area and magnitude of each negative clump which is, by its nature, a false peak detection. We can see in Fig. 6 that the inferred distribution of false clumps is found at small areas (peaking at ~1 kpc²) and faint luminosities (peaking at M_F814W ~ −17.4 mag) as expected.

Most of the negative detections are found at M_clump > −18 which is fainter than the peak for the fainter clumps that is observed in the magnitude distribution (see Fig. 10). This indicates that the smaller clumps with M_clump < −18 are mostly real physical clumps.

Fig. 6 Black points: distribution in areas and magnitudes for the negative clumps detected applying the same method on the inverted galaxy image. The contours depict the 68%, 90%, and 95% levels of the positive clump detections in our data. Most of the negative detections are found at Mclump > −18 which is fainter than the peak for the fainter clumps that is observed in the magnitude distribution (see Fig. 10). This indicates that the smaller clumps with Mclump < −18 are mostly real physical clumps.

3.6. Impact of the spatial resolution and morphological k-correction

We assumed that clumps are evenly detected across the rest-frame wavelength range that we probed at different redshifts. We tested this assumption by computing the same quantities on a version of the F814W images degraded to the same spatial resolution as the H-band CANDELS images, and on a sub-sample of galaxies for which H-band imaging from CANDELS is available. To degrade the I-band imaging data we re-scaled the image to match the pixel-scale of the H band and then convolved the image with a Gaussian kernel to match the PSF resolution of that image.

We find that the clumpy fraction in the degraded I-band images is found to be approximately three times smaller than what we find for the non degraded images. In the H-band imaging we find that the clumpy fractions are fully consistent with the values we get for the degraded I-band imaging data, possibly slightly lower. We therefore argue that the differences we find when comparing I- to H-band data can be explained mostly by the difference in spatial resolution which, as expected, merges clumps that are too close to be resolved. We note however that we find no evidence for the existence of a morphological k-correction that needs to be applied to f_clumpy.

3.7. Comparing with the Guo et al. (2015) and Shibuya et al. (2016) clump finding method

For comparison purposes, we have also implemented the same methodology as used in recent studies (e.g., Guo et al. 2015; Shibuya et al. 2016). As discussed at the start of this Sect. 3 the method used by Guo et al. (2015) and Shibuya et al. (2016) is based on selecting clumps with respect to the local image background. We implement the following algorithm:

• 1.

  SExtractor is run on a 7″ by 7″ image to obtain the segmentation map.

• 2.

  The segmentation map is dilated using a binary operation that effectively smooths and enlarges by a maximum of three pixels the outer edges of the segmentation map.

• 3.

  A small cutout of the original image is selected to adapt to the size of the segmentation map.

• 4.

  A boxcar filter (ten pixels width) is applied to the new cutout.

• 5.

  A residual image is produced by subtracting the filtered image to the original one.

• 6.

  The residual background flux is measured using a 3σ-clipping method.

• 7.

  All pixels that have fluxes smaller than 2σ of the calculated background are masked.

• 8.

  We run a new detection algorithm on the masked residual image to detect clumps.

All regions in the masked residual image that have more than five connected pixels are counted as single clumps. We also computed the magnitude of each clump using a three-pixel aperture.

To compute the final clump fraction with this method, we selected only galaxies for which GALFIT is converging into a reliable fit (see Ribeiro et al. 2016) and with an effective radius r_e > 0.2″ and axis ratio q > 0.5. We selected galaxies with stellar masses in an intermediate range (9.8 < log ₁₀(M_⋆/M_⊙) < 10.6). As in Guo et al. (2015) we also included single off-center clumps when computing the clumpy fraction. Finally, we only considered clumps which contain more than 8% of the galaxy light. We argue however that relying on GALFIT results for the final sample selection is limiting the potential for finding highly disturbed systems as well as very compact sources. It may also be that the structural parameters that we would rely on for that selection can be biased toward the brightest clump and as such not being representative of the entire galaxy, as we show in Ribeiro et al. (2016). Thus we believe that our methodology, which allows for a morphologically unbiased sample and includes a specific treatment of the apparent evolution caused by surface brightness dimming effects, can be complementary to that followed by Guo et al. (2015) and Shibuya et al. (2016) to provide a new perspective on the clumpy nature of high redshift galaxies. We compare the results obtained with this method to our results in Sect. 4.2.3.

Fig. 7 Number of galaxies as a function of the number of clumps detected within the selected isophote. We note that we include in this plot single clump galaxies which correspond to galaxies with no substructure or a single central clump. We plot this quantity for two redshift intervals each with its respective value of kp = 1.5 for 2 < z < 4.5 in gold (left panel) and kp = 5 for 4.5 < z < 6.0 in purple (right panel). On each panel we show the number of clumps for the VUDS targets.

4. The evolution of the clumpy fraction with redshift

4.1. The distribution of the number of clumps per galaxy

We compare in Fig. 7 the number of galaxies as a function of the number of detected clumps per galaxy in two different redshift ranges. We consider only clumps that reside inside the isophote belonging to the VUDS target when showing this histogram. We find that the fraction of multi-clump systems decreases as the number of clumps increases. Systems with only one clump (i.e., a single clump in a smooth galaxy at the observed resolution) are the dominant population (~55% at 2 < z < 4.5 and ~61% at 4.5 < z < 6.0). We find that the fraction of galaxies with only two clumps dominates the distribution of clumpy galaxies with ~25 ± 2.7% double clumps galaxies at 2 < z < 4.5 and ~21 ± 7.7% at 4.5 < z < 6.0. Galaxies with three clumps are about ~11 ± 1.6% at 2 < z < 4.5 and ~8 ± 4.1% at 4.5 < z < 6.0. Galaxies with four or more clumps amount to only ~7 ± 1.2% at 2 < z < 4.5 and ~9 ± 4.5% at 4.5 < z < 6.0. Several examples of galaxies with one, two, or more clumps may be found in Figs. 3−5, respectively.

4.2. The fraction of clumpy galaxies

We defined any source that has more than one clump detected within its isophote as a clumpy galaxy. To compute this fraction as a function of redshift, we considered only clumpy (and non-clumpy) galaxies within VUDS for which we have reliable information on the spectroscopic redshift. This ensures that we did not smear our results between redshift bins by including galaxies with uncertain or incorrect redshifts. We find that the stellar mass distribution of clumpy and non-clumpy galaxies are indistinguishable from one another, so that there is no particular population bias in defining clumpy galaxies as those with more than one clump. We stress that the clumpy fraction with our method depends on the selected isophote (value of k_p) and its dependency on a color correction, as allowed by the depth of available images.

We find that the clumpy fraction rises with redshift for all values of k_p that we have tested as well as when using near ultraviolet (NUV) as the base rest-frame band. Due to the asymmetric nature of most galaxies at high redshift (e.g., Huertas-Company et al. 2015; Curtis-Lake et al. 2016) and since deriving a center and a size can lead to different results when using different assumptions (e.g., Ribeiro et al. 2016) we did not include in our sample of clumpy galaxies those sources with an off-center clump as defined by, for example, Guo et al. (2015), Shibuya et al. (2016). The clumpy fraction evolution with redshift is shown in Fig. 8 and summarized in Table 2 and is discussed below.

4.2.1. Clumpy fraction in 2 < z < 4.5

At 2 < z < 4.5 we observe that the clumpy fraction remains in the range 35 to 55% at all redshifts. This fraction is observed to be rising from low to high redshifts. We also report the clumpy fraction for two different values of k_p showing that the absolute value does indeed depend on the isophote used. Increasing the value of k_p lowers the probability of detecting faint clumps and it is no surprise that the fraction of clumpy systems decreases as the value of k_p increases. We find, regardless of the value of k_p adopted, a factor ≳2 increase in the clumpy fraction from the lowest to the highest redshifts measured in our sample. The rising evolution with redshift of the clumpy fraction observed using our method is different from results obtained in other studies (e.g., Guo et al. 2015; Shibuya et al. 2016). A possible explanation for these differences is that the lower fraction of clumpy galaxies that we observe at z ~ 2−3 is the result of the evolving isophote threshold defined in Sect. 3, including luminosity evolution. This imposes a higher threshold at lower redshifts than studies using a fixed apparent luminosity threshold, resulting in a reduced area in which to search for clumps and therefore a lower number of clumps at the lower redshift end. Our choice of taking into account luminosity evolution in defining the search isophote is motivated by the objective of counting clumps within similar physical areas while stellar populations brighten in galaxies when moving from the highest redshifts in our sample down to the peak in star formation rate density. When selecting a fixed apparent isophotal level like in Guo et al. (2015) and Shibuya et al. (2016), the global brightening of galaxies enables to count clumps in larger galaxy areas making it more difficult to relate clump measurements at different redshifts.

Fig. 8 Fraction of clumpy galaxies within the VUDS sample as a function of spectroscopic redshift. Different large colored circles (gold, green, and purple) represent different values of the normalization kp. The smaller gray symbols represent values from the literature for clumpy galaxies: circles (Shibuya et al. 2015); squares (Guo et al. 2015); and for asymmetric galaxies: inverted triangles (Ravindranath et al. 2006); triangles (Conselice & Arnold 2009); octagons (Law et al. 2012); pentagons (Curtis-Lake et al. 2016). For kp = 5, the isophote at z < 4 is very high and as we only probe the brightest regions of galaxies we miss those clumps outside the search isophote, therefore lowering the fraction of clumpy galaxies. The empty purple circle at z ~ 5 is the extrapolation of the clumpy fraction as if it was measured at kp = 1.5 (see Sect. 4.2.2 for details).

4.2.2. Clumpy fraction 4.5 < z < 6

At z > 4.5 the depth of the data only allows the use of a higher value of k_p, as shown in Fig. 8. At these redshifts it is also necessary to use a brighter isophote limit for the clumps detection to avoid excessive contamination from noise. With these constraints we observe a clumpy fraction of 20−40% in the redshift range 4.5 <z < 6.

It is clear from our results at lower redshifts that if the depth of the data was enabling us to use a fainter isophote we would likely observe a higher clumpy fraction in this redshift bin too. Since we have computed the clumpy fraction at all redshifts for all values of k_p shown in Fig. 8 we attempted to compute a correction to estimate the clumpy fraction in 4.5 <z < 6 as if it was measured with k_p = 1.5. To do that we fitted a linear curve to the quantity $\frac{{\mathit{f}}_{\mathit{c,K}}}{{\mathit{f}}_{\mathit{c,}\mathrm{5}}}\mathrm{\left(}\mathit{z}\mathrm{\right)}$ which quantifies the underestimation of the clumpy fraction measured at k_p = 5 with respect to that measured at k_p = K. We then applied this correction factor measured at z = 5.25 to correct the clumpy fraction measured for k_p = 5.0. In doing so, we estimated that the clumpy fraction at 4.5 <z < 6.0 is ~64 ± 19% when corrected for k_p = 1.5.

4.2.3. Comparison of the clumpy fraction as derived from the Guo et al. (2015) method

We compared the results obtained with our method to results obtained using the Guo et al. (2015) method as defined in Sect. 3.7 (see Fig. 9). We find that for z < 3 we are in excellent agreement with the results from Guo et al. (2015) and Shibuya et al. (2016) when we used a method similar to theirs. At higher redshift we included our clumpy fraction measurement binned at 3 <z < 4 to minimize statistical errors. We find that our clumpy fraction is higher than what is reported by Shibuya et al. (2016) but it is compatible within ~2σ of our reported measurements. We thus conclude based on the bins for which we have larger samples hence more robust statistics that the main difference between our results and those reported by Guo et al. (2015) and Shibuya et al. (2016) is due to the different methods used to detect clumps in high redshift galaxies.

Fig. 9 Clumpy fraction in the VUDS sample as computed using the method described by Guo et al. (2015; red squares) and the results by Guo et al. (2015) as blue diamonds and Shibuya et al. (2016) as gray circles which use the same methodology.

5. Physical properties of clumps

5.1. The relation between clump luminosity and area

As stated in Sect. 3.3, we derived areas and absolute magnitudes for each clump by computing the flux and area within each group of pixels associated with the local maximum that defines the clump. We show in Figs. 10 and 11 the tight correlation that we observe between clump magnitude and clump area. This is expected by construction since the larger the area (i.e., the number of pixels of the clump) the more likely it is for it to be brighter. In this figure we exclude galaxies with a single clump detection.

We can observe in Fig. 10 that the distribution of clump areas is quite broad, and that the magnitude distribution has an asymmetric shape. We fitted each of these distributions with a two Gaussian model with six free parameters: two amplitudes, two centers and two widths, one for each Gaussian curve. The fit was done using the Levenberg-Marquardt minimization algorithm implemented in scipy. For clumps at 2 <z < 4.5 we observe two peaks in the magnitude distribution at M_NUV = −18.1 ± 0.1 and M_NUV = −19.7 ± 0.1. In the area distribution the two peaks are found at log (T) = 0.25 ± 0.05 and log (T) = 0.78 ± 0.03, which correspond to ~1.8 and ~6.0 kpc². We test the statistical robustness of a possible bimodality using the Hartigan’s dip test (Hartigan & Hartigan 1985). We find that the unimodal distribution hypothesis is highly favored (at ~80% level). We do find a deviation from a normal distribution computed from the Shapiro & Wilk (1965) test that rejects the normal distribution hypothesis at a 99% level. At 4.5 <z < 6.0 the low number of clumps does not allow us to test for a bimodal distribution (see Fig. 11).

We conclude from this analysis that a broad range of clump luminosities and sizes is observed. Clumps are confidently identified with areas as small as 1 kpc² and faint luminosities M_clump,FUV ~ −17, or with sizes as large as 15 kpc² and bright luminosities M_clump,FUV ~ −21. This diversity in clump properties may be linked to different processes of clump formation producing bright and large clumps preferentially with M_clump,FUV < −18.5 and area T_clump > 4 kpc², or faint and small clumps with M_clump,FUV > −18.5 and area T_clump < 4 kpc², producing the continuous range of observed properties. This is discussed in Sect. 6.

Fig. 10 Absolute magnitude in the rest-frame FUV and area distribution for clumps found in galaxies with Nc ≥ 2 detected at 2 <z < 4.5 with kp = 1.5. Each clump is color coded by the number of clumps found in its parent galaxy. The color scheme is found at the bottom right of the main panel. The black solid line in each of the auxiliary panels shows a fit of a two Gaussian model to the observed distribution of areas and magnitudes of the clumps. We observe that it is more probable to find the brightest and largest clumps in two-clump systems. The contours show the smoothed distribution of the fraction of two-clump systems at levels of 10%, 25%, and 50%. It is clear that the majority of two-clump galaxies (with respect to all clumpy systems) is found in brighter and larger clumps.

5.2. Bright and large clumps are in two-clump galaxies

A key piece of information shown in Fig. 10 is the number of clumps in the galaxy that each clump belongs to. The larger and brightest clumps are most likely to be found in two-clump galaxies and as we move toward fainter and smaller clumps we observe that they are more commonly found in galaxies with a larger number of clumps. This information is summarized in the contour lines of Fig. 10 showing the distribution of the fraction of clumps in two-clump galaxies as a function of clump area and absolute magnitude. More than 50% of the brightest and largest clumps are found in two-clump galaxies and the probability of finding two clumps in a galaxy decreases for fainter clump magnitudes and smaller areas (down to ~10%). The information on how the fraction of clumps in two-clump galaxies evolves with redshift is found in Table 3.

5.3. Luminosity ratios of two-clump galaxies

To classify each two-clump galaxy we compute the luminosity ratio between the two clumps, L_f/L_b, with L_f being the luminosity of the fainter clump and L_b that of the brighter of the two. The resulting fractions as a function of redshift can be found in Table 3. The fraction of two-clump galaxies is classified in major (1/4 <L_f/L_b < 1) and minor (L_f/L_b < 1/4) systems, for galaxies which have at least one clump with M_clump < −19.1. Statistical uncertainties become large because of the reduced number of galaxies for which we are able to compute this quantity.

5.4. Clump stellar masses

We compute the stellar mass of clumps with a simple assumption that the mass to light ratio of clumps is the same as its host galaxy. This hypothesis allows us to estimate the stellar mass of each clump by computing the clump to galaxy luminosity ratio and applying this ratio to compute the stellar mass of a clump from the total stellar mass of the galaxy. Clearly, this estimate is a first order approximation as stellar clumps are likely to be on average more star-forming than the galaxy as a whole if they are strongly star-forming regions, and we might therefore underestimate the stellar masses of clumps.

In Fig. 12 we compare the distribution of stellar masses in clumps for three different samples. A first sample is considering clumps that are in clump systems with at least three clumps (hereafter multi-clump systems). Another stellar mass distribution is for the sample of clumps in two-clump galaxies. These two samples are compared to the distribution of stellar masses in single clumps galaxies. We note that for single clump systems, a clump can correspond the entire galaxy, explaining why we find clumps with stellar masses up to log ₁₀(M_⋆/M_⊙) ~ 11.

From Fig. 12 it is striking to see that the clumps in two-clump galaxies have stellar masses that peak at higher values (around log ₁₀(M_⋆/M_⊙) ~ 9.4) than clumps in multi-clump systems. The stellar masses of these clumps appear to be larger than what is expected for giant clumps identified from galaxy disk simulations or observations (see e.g., Elmegreen et al. 2009; Genzel et al. 2011; Guo et al. 2012; Bournaud 2016), with such clumps having lower stellar masses in the range 8 ≲ log ₁₀(M_⋆/M_⊙) ≲ 9. On the other hand, the bulk (~55%) of the stellar mass distribution in multi-clump-systems is consistent with this lower mass range. We have tested the similarity between the distributions of stellar masses in two-clump galaxies and in multi-clump systems using a Kolmogorov-Smirnov test and the null hypothesis (that both distributions are the same) is discarded at the 95% confidence level. This difference persists when considering higher values of k_p. We conclude that these two populations are genuinely different, and discuss the possible reasons for these differences in Sect. 6.

Fig. 11 Absolute magnitude in rest-frame FUV and area distribution for clumps found in galaxies with Nc ≥ 2 detected at 4.5 <z < 6.0 with kp = 5. Each clump is color coded by the number of clumps found in its parent galaxy. The color scheme is found at the bottom right of the main panel. The contours shows the distribution of fraction of two-clump systems at levels of 10% and 25%.

Fig. 12 Clump stellar mass distributions for three galaxy subsamples. Black, gold, and filled pink histograms respectively indicate the distribution for galaxies with one clump, two clumps, and three or more clumps.

6. Discussion: a merger origin for the bright and large clumps in two-clump galaxies

Our study brings to light two important morphological properties of high redshift galaxies and their evolution from z ~ 6 to z ~ 2: the fraction of clumpy galaxies as well as the distribution of the number of clumps per galaxy. Both of these properties are a direct consequence of the evolutionary path followed by galaxies, and the physical processes driving this evolution at these early cosmic times.

We find that the fraction of clumpy galaxies in the star-forming population rises from ~ 34 ± 8% at z ~ 2.25 to 54 ± 9% at z ~ 4 for galaxies selected in the UV rest-frame with stellar masses log ₁₀(M_⋆/M_⊙) > −0.204 × (z−4.5) + 9.35, and identifying clumps brighter than M_clump,NUV < −18 and larger than ~2 kpc². Our method has the benefit of producing counts that can be readily compared over a large redshift range, revealing that the clumpy fraction remains high over the redshift range 2 <z < 6.

These results on the evolution of the clumpy fraction are apparently at odds with the decreasing fraction of clumpy galaxies reported by Guo et al. (2015) and Shibuya et al. (2016), but can be easily understood from the difference in the surface brightness limits used to identify clumps. The surface brightness limit evolving with redshift that we use seems to be the key difference with Guo et al. (2015) and Shibuya et al. (2016). We made the choice to use an isophote evolving with redshift (Eq. (2)), as motivated by the general luminosity evolution with redshift. In comparison, Guo et al. (2015) and Shibuya et al. (2016) use a fixed isophote level, independent of redshift. This implies that we use smaller search areas at the lowest redshifts (as isophote levels are higher in our method) and higher search areas at the highest redshifts (lower isophote levels in our method). This implies in turn that we count less clumps at low redshifts, and more clumps at high redshifts.

Other more subtle differences in the definition of a clumpy galaxy and in the sample used to compute the clumpy fraction may further contribute to the differences observed between our study and others. We define a clumpy galaxy as a galaxy with two or more detected clumps within its limiting isophote, whereas Guo et al. (2015) and Shibuya et al. (2016) define it when at least one off-center clump is identified, and we should then expect a lower clumpy fraction from our measurements. Another important difference resides in our sample selection which does not apriori exclude galaxies because of their morphology (Guo et al. 2015; Shibuya et al. 2016, impose an effective radius r_e > 0.2″ and axis ratio q > 0.5 cut). As we included in our sample small compact sources, this might further increase the number of non-clumpy galaxies and decrease the clumpy fraction, while including galaxies with low values of q like elongated clump chain galaxies (see e.g., Elmegreen et al. 2004, 2008) would rather increase it.

We therefore argue that the overall lowering of the clump fraction with increasing redshift as reported by previous studies may simply be due to the bias produced by using the same observed surface brightness at all redshifts when identifying clumps: because of the limits in detection at the highest redshifts one probes only the brighter and rarer clumps, while at the same observed limiting surface brightness one finds more numerous fainter clumps, hence mimicking a clump fraction decrease with redshift. This important methodological difference leads to results significantly different from previous literature results, and we believe that our method improves on earlier attempts with a more physically motivated method. Our report of a constant or slightly increasing clump fraction over 2 <z < 6 is consistent with the large range of total sizes of galaxies that we report in a previous work (Ribeiro et al. 2016), and indicates that the morphology of galaxies remains varied and complex from the end of reionization to the peak in star formation activity.

The other particularly important result from our paper is from the analysis of the distribution of the number of clumps per galaxy, were we observe a large fraction (~21−25%) of galaxies with only two clumps. The clumps in this fraction of the galaxies are the brightest and largest as shown in Figs. 4, 7, 10, and 11. We can identify only a few possibilities to form the clumps in these two-clump galaxies: they could be either (1) the result of disk fragmentation and VDI; (2) two individual galaxies in the final merging process of a galaxy pair; or (3) the chance projection of a galaxy on the line of sight, physically unrelated. We note that we are able to detect clumps that are one to two magnitudes fainter than what we find for bright clumps in two-clump systems. This implies that if additional clumps were to be found with deeper imaging they would be much fainter than what we find and thus we would still observe a large fraction of two-clump systems detected above a given luminosity limit.

Galaxies with VDI-induced clumps have been studied in detailed hydrodynamical simulations which show how fragmentation in a proto-galactic disk produces several clumps, visible over timescales ≲1 Gyr, before migrating and merging in a galaxy central core (e.g., Bournaud et al. 2007, 2014; Dekel et al. 2009; Genel et al. 2012; Oklopčić et al. 2017). It seems unlikely that a physical scenario following disk fragmentation and VDI would preferentially produce only two bright clumps (F. Bournaud, priv. comm.). Cloud fragmentation due to disk instabilities are more likely to give rise to a higher number of clumps, with the number of clumps reducing over the dynamical timescale of clump merging, as shown in simulations (Mandelker et al. 2014). If VDI is the only process capable of producing clumps one would rather expect to observe more galaxies with three or more clumps than two-clump galaxies, contrary to what we find (Fig. 6). It is however possible that some of the two-clump galaxies that we observe are the end result of clumps formed through VDI after migration and merging of smaller clumps. The typical timescale for clump migration to a galaxy center is about 500 Myr (e.g., Förster Schreiber et al. 2011; Guo et al. 2012; Bournaud 2016) which is typically one fourth of the redshift interval we probe from 2 <z < 4.5. However, the period on which we would observe two resolved clumps should be even shorter, further reducing the fraction of large and bright clumps if produced by VDI. Additionally, we would need to take into account the survival rate of large clumps for which some simulations suggest even shorter lifetimes (e.g., Genel et al. 2012). We therefore infer that the large fraction of two-clump galaxies with large and bright clumps cannot be produced by VDI processes alone.

Another interpretation is that these two-clump galaxies are direct evidence for on-going merger events. Major or minor mergers do indeed produce a bi-modal light distribution with each of the merging galaxy counted as a clump. Merging is observed to be a significant process with a high fraction of ~15−20% of galaxies involved in major mergers at redshifts z ~ 1 and up to z ~ 4 (e.g., Le Fèvre et al. 2000; de Ravel et al. 2009; López-Sanjuan et al. 2013; Tasca et al. 2014a). The bright clumps in two-clump galaxies and a luminosity ratio L_f/L_b > 1/4 have absolute magnitudes −19 >M_NUV > −22, producing together a total luminosity equivalent to that of the median of the star-forming population we are probing. Two-clump galaxies with one bright clump and L_f/L_b < 1/4 may include minor merger events participating to the build-up of stellar mass along cosmic time.

If we assume that clumps which have L_f/L_b > 1/4 are representative of major galaxy mergers we find that 24−40% of two-clump galaxies are then classified as major merger events, remaining roughly constant at z ≳ 3. As two-clump galaxies represent ~25% of the total population, using the fraction of bright two-clump galaxies (Table 3) we find a total merger fraction slightly decreasing with redshift from ~18% at z ~ 3 to ~10% at z ~ 5. Using those systems with L_f/L_b > 1/4 we find a major merger fraction in the range 6.5 to 11%. Similarly, using the fraction of two-clump galaxies with L_f/L_b < 1/4 we find a minor merger fraction of 2.7 to 10.5%.

The fraction of merging systems derived above is obtained when considering close merging pairs alone, hence late stage mergers, as the median clump separation is ~3 kpc (~0.4′′). To expand this to a total merger fraction we need to account for the fraction of galaxy pairs within a 20 kpc search radius, a distance for which the probability of galaxies to merge in a short timescale is expected to be high (e.g., Kitzbichler & White 2008). From a pair count corrected for projections along the line of sight (Sect. 3.4), we find a rising pair fraction of ~15% at z ~ 2.3 reaching ~38% at z ~ 5. One interesting trend is that while the pair fraction is increasing at higher redshifts, the opposite trend is observed for the merger fraction derived from clumps statistics (f_2c in Table 3). One may infer that we are observing an evolving merger rate measured in two stages, the first and more common at the higher redshifts is when galaxies have yet to merge, and the second, as observed from close bright clumps counts, representing the final phases of merging before coalescence, which is more common at lower redshifts. When combining these two merging indicators we find that the major merger rate remains roughly constant at ~20% over the redshift range 2 <z < 6.

There is a low probability that our observations of two-clump galaxies might be consequence of line of sight contamination from lower redshift galaxies which might then be confused with high redshift clumps. With a small median separation of ~0.4′′, the probability that these two-clump galaxies with at least one bright M_clump,NUV < −19.1 are the result of random line of sight contamination separations is very small, at less than 3% (see Sect. 3.4 for details), even when dealing with faint objects with i_AB ~ 27−28. We then conclude that most of these two-clump galaxies are real physical systems at the same redshift.

As a conclusion, a straightforward interpretation of galaxies with two bright clumps is that this population predominantly includes galaxies caught in the act of a merging event. This interpretation is further reinforced when looking at the typical stellar masses of clumps. We find that for those two-clump systems with at least one bright component, the majority (~77%) have stellar masses greater than log (M_⋆/M_⊙) > 9 a limit beyond which clumps induced by VDI are not observed in simulations (see e.g., Elmegreen et al. 2009; Guo et al. 2012; Tamburello et al. 2015; Bournaud 2016). On the other hand, ~55% of clumps in multi-clump systems have stellar masses lower than log (M_⋆/M_⊙) < 9, which might likely result from VDI as well as minor mergers. We find that the typical masses of the clumps in two-clump systems are more alike the masses of galaxies observed in major merging pairs with larger separations (López-Sanjuan et al. 2013; Tasca et al. 2014a).

From the analysis of clump properties presented in this paper, we therefore conclude that we are likely witnessing two different modes of clump formation happening in parallel at the same cosmic times. Small and faint clumps are likely the result of both VDI clumps production and minor mergers. At the same cosmic time, we do observe an important fraction (in ~20−25% of all galaxies) of two-clump galaxies with large and bright clumps which are most likely major mergers.

While our conclusions are based on indirect evidence, direct observational evidence will be required to confirm the physical origin of the two populations of clumps that we have identified. The definite proof that the bright two-clump galaxies are indeed major merging systems would require dynamical information on the velocity field, with the relative velocity separation of the two bright clumps. Some of these merging systems are already confirmed from spectroscopic evidence in the VUDS survey (Tasca et al. 2014a), and more evidence of a high merger fraction will be presented in a forthcoming paper (Le Fèvre et al., in prep.). Likewise, velocity maps of multiple clump systems could bring some more evidence in support of clump formation via instabilities, although the dynamical signature of a VDI process may be rather subtle (Bournaud et al. 2009) and may escape current or even upcoming observational capabilities. Only the result of the VDI process in the form of a clump or several clumps, as shown in simulations, is accessible to observations today as reported here. Identifying minor mergers within the population of small and faint clumps is also important but again will require significantly improved observational capabilities. Further work is therefore needed to consolidate the picture we propose, but our results point out the diversity of clump properties, which signals a likely diversity in their formation process.

7. Summary

The results presented in this paper can be summarized as follows:

• The search for clumps in distant galaxies images criticallydepends on the depth of the images defining the search area, onhow clumps are defined, whether clumps must be embeddedwithin a common isophote, and on the selection of the galaxypopulation. We used a new algorithm for detecting clumps, inwhich clumps are defined as groups of connected pixelsassociated with the same local maximum within a galaxyisophote evolving with redshift following surface brightnessdimming and luminosity evolution, and corrected for colordifferences between the different rest-frame bands observed as afunction of redshift.

• We computed the number of clumps per galaxy in a sample of 1242 galaxies with 2 < z < 4.5 and 96 galaxies with 4.5 < z < 6 with confirmed spectroscopic redshifts from the VIMOS Ultra Deep Survey and in the COSMOS and ECDFS areas where deep HST imaging in the F814W filter is available from the COSMOS and CANDELS surveys.

• The fraction of clumpy galaxies, with two or more clumps and corrected for line of sight chance projections, is found to remain in the range 35 to 55% over 2 < z < 6, possibly increasing with redshift.

• We find that at 2 <z < 4.5 and excluding single clump galaxies, the dominant galaxy population has two clumps (~25% of all galaxies) followed by galaxies with three, and more than three clumps (~11, ~7% respectively). The same trend is found at higher redshifts (4.5 <z < 6) where galaxies with two clumps dominate (~21%) followed by galaxies with three, and more than three clumps (~8%, ~9% respectively).

• We find a large range of clump properties, with a continuous, possibly bi-modal, distribution of the luminosity and area of clumps with a population of large and bright clumps distinct from a population of small and less luminous clumps. The population of large and bright clumps has typical absolute magnitudes M_FUV ~ −19.7 and areas T ~ 6.0 kpc². In contrast, the population of small and less luminous clumps has typical absolute magnitudes M_FUV ~ −18.1 and areas T ~ 1.8 kpc². We associate these two populations to the existence of two physical processes, mergers and disk instabilities, that explain the observed properties of clumps.

• The faint and small clumps are preferentially found in multiple clump galaxies with N_clump ≥ 3. They are comparable in luminosity and stellar mass (10⁸−10⁹ M_⊙) to the large star clusters seen forming in numerical simulations from violent disk instabilities, consistent with findings from other studies (e.g., Elmegreen et al. 2009; Guo et al. 2012; Bournaud 2016). It is also likely that a fraction of these clumps is related to minor merging events.

• The bright and large clumps are preferentially found in two-clump galaxies. These bright clumps have luminosities and stellar masses ~10⁹−10¹⁰ M_⊙, larger than expected from VDI processes. It appears as unlikely that, if clumps were solely produced by VDI, two-clump galaxies would strongly dominate the population of clumpy galaxies. We rather conclude that the bright and large clumps in two-clump systems have properties similar to those of galaxies in merging pairs (López-Sanjuan et al. 2013; Tasca et al. 2014a). Using the fraction of galaxies with bright clumps we infer a major merger fraction going from 18% at z ~ 3 to 10% at z ~ 5. When combined with the pair fraction of galaxies as measured in projection within a 20 kpc radius, we find a major merger rate of ~20% comparable to other estimates at similar redshifts (López-Sanjuan et al. 2013; Tasca et al. 2014a).

In quantifying the properties of clumps in star-forming galaxies with 2 < z < 6 we therefore conclude that we are witnessing the effects of two different galaxy assembly processes working in parallel at a cosmic time of major galaxy assembly. Clumps are likely to be the result of both violent disk instabilities with an in-situ clump formation mode, as well as from major and minor merging of galaxies assembling matter coming from different origins. On the one hand gas accretion is expected to feed the buildup of disks which then become unstable and produce clumps that will build galaxy bulges (e.g., Dekel et al. 2009; Bournaud et al. 2009; Dekel & Krumholz 2013; Bournaud 2016), a process which remains to be confirmed by further observational evidence at these high redshifts (e.g., Bouché et al. 2013). The population of large and bright merging clumps formed ex situ from the galaxy they are assembling into confirm mounting evidence for the importance of major mergers in the early build-up of galaxies, events which contribute the most to the stellar mass growth of a galaxy (López-Sanjuan et al. 2013; Tasca et al. 2014a; Mandelker et al. 2014).

Both disk formation processes and early merging events are therefore important to understand both the evolution of the cosmic star formation rate and the build-up of stellar mass of galaxies throughout cosmic time, and explaining key observables such as the evolution of the specific star formation rate, and its possible flattening at redshifts z > 2 (e.g., Dekel & Mandelker 2014; Tasca et al. 2015; Faisst et al. 2016; Mármol-Queraltó et al. 2016). We point out that in this study we did not discuss the properties of galaxies with a single bright clump, which dominate the population of star-forming galaxies. This population will be the subject of forthcoming papers.

Testing the evolution of the number of clumps from larger simulations of galaxies in complete cosmologically representative volumes would bring further insight into the relative contribution of these different physical processes in building-up galaxies along cosmic time. On the observational side, since the clumpy fraction and the number of clumps identified strongly depend on the spatial resolution and depth of images, and as galaxies with large sizes are observed up to large redshifts (Ribeiro et al. 2016), it is likely that future facilities like JWST and extremely large telescopes with adaptive optics may well identify more clumps. Spatially resolved information from multi-band high resolution observations and/or IFU data will be necessary to compute velocity fields and stellar ages of the large and bright clumps and confirm or not a scenario where these are witnesses of the final stages of merging pairs of galaxies. The high fraction of bright two-clump galaxies and the wide range of clump properties reported here must be taken into account when planning for future surveys aimed at a complete census of the properties of the high redshift population.

Acknowledgments

We thank the anonymous referee for his/her comments and suggestions which helped improve this manuscripts. This work is supported by funding from the European Research Council Advanced Grant ERC–2010–AdG–268107–EARLY and by INAF Grants PRIN 2010, PRIN 2012 and PICS 2013. R.A. acknowledges support from the ERC Advanced Grant 695671 “QUENCH”. This work is based on data products made available at the CESAM data center, Laboratoire d’Astrophysique de Marseille. B.R. gratefully thanks Ana Afonso for late-night and helpful discussions, style tips and English typo corrections. This research made use of Astropy, a community-developed core Python package for Astronomy (Astropy Collaboration et al. 2013); Numpy & Scipy (Walt et al. 2011); Matplotlib (Hunter 2007).

References

All Tables

All Figures

	Fig. 1 Stellar masses and star formation rates for the galaxies in the selected sample. In each panel, the light gray points show all VUDS galaxies with measured redshift in each bin. The colored circles show the 1242 galaxies in the selected stellar mass range and with correct redshift probability of >75% (see text for more details).
In the text

Fig. 2 Exemplification of the method for detecting galaxy clumps. On the left we show the number of disconnected regions with more than ten pixels as a function of flux level cut. When f = 1 the isophote limit represents the limiting area of the galaxy. On the right we show the contours corresponding to three different flux level cuts at 20%, 50% and 100% of total flux that are identified in the same color as the vertical lines in the left panel.

In the text

	Fig. 3 Twelve examples of galaxies for which only one clump was detected. The white circle shows the 20 kpc projected search radius. The contours represent the detection isophote of each galaxy identified within this search radius. Each circle shows the position of a detected clump within the isophote of the same color. VUDS targets are identified with red contours and circles.
In the text

Fig. 4 Twelve examples of galaxies with two clumps detected. The top six examples show cases for which at least one bright clump was detected whereas the bottom six examples show pairs of faint clumps. The white circle shows the 20 kpc projected search radius. The contours represent the detection isophote of each galaxy identified within this search radius. Each circle shows the position of a detected clump within the isophote of the same color. VUDS targets are identified with red contours and circles.

In the text

	Fig. 5 Twelve examples of galaxies with more than two clumps detected. The white circle shows the 20 kpc projected search radius. The contours represent the detection isophote of each galaxy identified within this search radius. Each circle shows the position of a detected clump within the isophote of the same color. VUDS targets are identified with red contours and circles.
In the text

Fig. 6 Black points: distribution in areas and magnitudes for the negative clumps detected applying the same method on the inverted galaxy image. The contours depict the 68%, 90%, and 95% levels of the positive clump detections in our data. Most of the negative detections are found at Mclump > −18 which is fainter than the peak for the fainter clumps that is observed in the magnitude distribution (see Fig. 10). This indicates that the smaller clumps with Mclump < −18 are mostly real physical clumps.

In the text

Fig. 7 Number of galaxies as a function of the number of clumps detected within the selected isophote. We note that we include in this plot single clump galaxies which correspond to galaxies with no substructure or a single central clump. We plot this quantity for two redshift intervals each with its respective value of kp = 1.5 for 2 < z < 4.5 in gold (left panel) and kp = 5 for 4.5 < z < 6.0 in purple (right panel). On each panel we show the number of clumps for the VUDS targets.

In the text

Fig. 8 Fraction of clumpy galaxies within the VUDS sample as a function of spectroscopic redshift. Different large colored circles (gold, green, and purple) represent different values of the normalization kp. The smaller gray symbols represent values from the literature for clumpy galaxies: circles (Shibuya et al. 2015); squares (Guo et al. 2015); and for asymmetric galaxies: inverted triangles (Ravindranath et al. 2006); triangles (Conselice & Arnold 2009); octagons (Law et al. 2012); pentagons (Curtis-Lake et al. 2016). For kp = 5, the isophote at z < 4 is very high and as we only probe the brightest regions of galaxies we miss those clumps outside the search isophote, therefore lowering the fraction of clumpy galaxies. The empty purple circle at z ~ 5 is the extrapolation of the clumpy fraction as if it was measured at kp = 1.5 (see Sect. 4.2.2 for details).

In the text

	Fig. 9 Clumpy fraction in the VUDS sample as computed using the method described by Guo et al. (2015; red squares) and the results by Guo et al. (2015) as blue diamonds and Shibuya et al. (2016) as gray circles which use the same methodology.
In the text

Fig. 10 Absolute magnitude in the rest-frame FUV and area distribution for clumps found in galaxies with Nc ≥ 2 detected at 2 <z < 4.5 with kp = 1.5. Each clump is color coded by the number of clumps found in its parent galaxy. The color scheme is found at the bottom right of the main panel. The black solid line in each of the auxiliary panels shows a fit of a two Gaussian model to the observed distribution of areas and magnitudes of the clumps. We observe that it is more probable to find the brightest and largest clumps in two-clump systems. The contours show the smoothed distribution of the fraction of two-clump systems at levels of 10%, 25%, and 50%. It is clear that the majority of two-clump galaxies (with respect to all clumpy systems) is found in brighter and larger clumps.

In the text

	Fig. 11 Absolute magnitude in rest-frame FUV and area distribution for clumps found in galaxies with Nc ≥ 2 detected at 4.5 <z < 6.0 with kp = 5. Each clump is color coded by the number of clumps found in its parent galaxy. The color scheme is found at the bottom right of the main panel. The contours shows the distribution of fraction of two-clump systems at levels of 10% and 25%.
In the text

	Fig. 12 Clump stellar mass distributions for three galaxy subsamples. Black, gold, and filled pink histograms respectively indicate the distribution for galaxies with one clump, two clumps, and three or more clumps.
In the text