Free Access
Issue
A&A
Volume 659, March 2022
Article Number A16
Number of page(s) 31
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202141651
Published online 25 February 2022

© ESO 2022

1. Introduction

The reionization of the Universe is an outstanding problem that remains unresolved. While the time scales over which reionization ended are well established around redshift z ∼ 6 (e.g., Mason et al. 2019; Yang et al. 2020; Paoletti et al. 2020), the dominant sources of photons responsible for the transformation of the dominant neutral hydrogen into a mostly ionized medium have yet to be determined. Faint low-mass star-forming galaxies are considered the candidates leading reionization in this era due to their large number density and weak gravitational potential, favouring the strong and effective feedback needed to open low HI density paths for photons to escape (e.g., Wise et al. 2014; Robertson et al. 2015; Bouwens et al. 2016; Finkelstein et al. 2019). However, the contribution of additional sources with higher ionizing photon efficiency, such as luminous, massive starburst galaxies (e.g., Naidu et al. 2020; Endsley et al. 2021) and active galactic nuclei (AGN, e.g., Grazian et al. 2018) might have a significant contribution (e.g., Dayal et al. 2020, and references therein).

A detailed characterization of the rest-frame ultraviolet (UV) spectra of star-forming galaxies (SFG) at z ≳ 6 is essential for understanding their ionization properties and, thus, shedding new light onto the reionization process. In the past few years, deep near-infrared (NIR) observations of galaxies during the reionization era have reported the presence of UV emission lines with unusually high equivalent widths (EW), such as CIVλ1550 Å, HeIIλ1640 Å, OIII]λ1663 Å, and CIII]λ1908 Å (see Stark 2016, for a review).

In particular, UV nebular lines encode precious information on the physical conditions of the ionized gas in galaxies. Different photoionization models that are used to understand the role of age, ionization parameter, metallicity, and dust on the emergent UV nebular lines struggle to explain their origin and strength (e.g., Jaskot & Ravindranath 2016; Gutkin et al. 2016; Feltre et al. 2016; Nakajima et al. 2018a; Byler et al. 2018; Hirschmann et al. 2019). However, thus far, all models require the presence of hard radiation fields that are capable of reproducing the observed UV emission lines with high ionization potentials, which also leads to more extreme ionization conditions in the interstellar medium (ISM; e.g., high EWs and line ratios). Constraining the available models with large and representative samples of emission line galaxies is therefore needed to improve our understanding of the physical mechanisms producing UV emission lines, thus paving the way for future extensive studies of galaxies at z > 6.

Emission lines are relevant not only to the understanding of the physical conditions governing these early galaxies but they also provide a tool for their spectroscopic redshift identification. This is especially relevant at z > 6, where absorption features are weak and a significant drop in the number of galaxies with Lyα emission, often the strongest emission line in the UV, is observed due to the sharp increase of absorption by a predominantly neutral intergalactic medium (IGM; e.g., Fan et al. 2006; Pentericci et al. 2014; Cassata et al. 2015; Fuller et al. 2020).

One of the best alternatives to Lyα is the CIII] doublet – a combination of [CIII]λ1906.68 Å, a forbidden magnetic quadrupole transition and CIII]λ1908.73 Å, a semiforbidden electro-dipole transition (here, we are referring to vacuum wavelength), with an ionizing potential of 24.4 eV. This doublet (hereafter cited as CIII]λ1908 or CIII], for simplicity) is typically the brightest UV metal line in star-forming galaxies at intermediate redshift (z ∼ 3, e.g., Shapley et al. 2003) and has been proposed as an alternative in the search for galaxies in the reionization era (Stark et al. 2014). Some searches have been successful, reporting high EW(CIII]) in the observed spectra of galaxies at z > 6 (Stark et al. 2015, 2017; Mainali et al. 2017; Laporte et al. 2017; Hutchison et al. 2019), but other studies have failed to detect the line in intensively star-forming systems (Sobral et al. 2015; Schmidt et al. 2016).

At low-z, studying UV emission lines requires space-based spectroscopy. Studies using Hubble Space Telescope (HST) observations of relatively small samples have shown that strong CIII] emission is generally present in the spectra of local low-metallicity galaxies (e.g., Garnett et al. 1995; Leitherer et al. 2011; Berg et al. 2016, 2018; Senchyna et al. 2017; Ravindranath et al. 2020). However, the characterization of larger samples spanning a wider range of properties, such as stellar masses, star formation rates (SFR), and metallicities, requires a stronger observational effort that has precluded studies with statistical significance. This is different at 2 < z < 4, where both Lyα and CIII] are redshifted into the optical and can be probed over larger samples with ground-based 8−10 m-class telescopes. Also, galaxies at z > 2 are likely to be more similar to those at z > 6 (see Shapley 2011, for a review). Several studies now routinely report CIII] emission (along with other strong emission and absorption lines) in galaxies at cosmic noon, either from small samples of relatively bright galaxies (Steidel et al. 1996, 2014; Erb et al. 2010; Amorín et al. 2017; Du et al. 2020), fainter gravitationally lensed galaxies (Pettini et al. 2000; Christensen et al. 2012; Stark et al. 2014; Rigby et al. 2015; Berg et al. 2018; Vanzella et al. 2021), or high signal-to-noise (S/N) stacks from larger galaxy samples (Steidel et al. 2001; Shapley et al. 2003; Le Fèvre et al. 2019; Nakajima et al. 2018b; Feltre et al. 2020).

Deep surveys such as the VIMOS Ultra Deep Survey (VUDS, Le Fèvre et al. 2015) or the MUSE Hubble Ultra Deep Field Survey (HDFS, Bacon et al. 2017) have recently studied large samples of CIII] emitters. Le Fèvre et al. (2019) showed that only 24% of the VUDS galaxies at 2.4 < z < 3.5 show CIII] emission; however, it is only in ∼1% that this emission is as intense as the values found at z > 6 and that is: EW(CIII]) ≳ 10−20 Å. Amorín et al. (2017) showed that extreme CIII] emitters at 2 < z < 4 in VUDS are very strong Lya emitters (LAEs), characterized by very blue UV spectra with weak absorption features and bright nebular emission lines. These galaxies present high excitation, low metallicities, and low carbon-to-oxygen (C/O) abundance ratios, similar to the common values expected in most of the galaxies during the first 500 Myr of cosmic time.

Based on the stacking of a large sample of LAEs from the MUSE HDFS, Feltre et al. (2020) found that the mean spectra of LAEs with larger Lyα EW, fainter UV magnitudes, bluer UV spectral slopes, and lower stellar masses show the strongest nebular emission. Maseda et al. (2017) arrived at similar conclusions for a sample of 17 CIII] emitters at z ∼ 2 in the MUSE HDFS. For these galaxies, they found a correlation between EW(CIII]) and EW([OIII]λ5007), linking the properties of the stronger CIII] emitters to those of the so-called Extreme Emission Line Galaxies (EELGs, e.g., Maseda et al. 2014; Amorín et al. 2015). These are low-metallicity starbursts defined by their unusually high EW([OIII]λ5007) ≳ 200 Å (see also, Tang et al. 2021; Matthee et al. 2021). At lower z, the Green Pea galaxies (Cardamone et al. 2009; Amorín et al. 2010, 2012) are EELGs and they include objects for which Lyman continuum leakage has been directly measured (Izotov et al. 2016; Guseva et al. 2020; Wang et al. 2021). These galaxies show prominent UV nebular lines, including high EW CIII] (Schaerer et al. 2018; Ravindranath et al. 2020).

While EELGs are likely analogs of the bright-end of reionization galaxies, the more common population of normal, main-sequence SFGs showing moderate or low EW(CIII]) still needs to be fully characterized at z > 2. This requires large and very deep samples achieving sufficiently high S/N spectra for detecting and studying the fainter CIII]-emitters. The present work is the first of two papers aimed at exploiting the unprecedented ultra deep spectra provided by the VANDELS survey (McLure et al. 2018; Pentericci et al. 2018) to assemble a large unbiased sample of main-sequence star-forming galaxies CIII] emitters at 2 < z < 4 and to characterize their main physical properties as a function of their UV line emission and chemical abundances.

The metal content of galaxies not only has a crucial role in the production and strength of nebular lines, but it is also sensitive to their star formation activity and to the presence of outflows, gas stripping, and dilution resulting from inflow of pristine gas (Maiolino & Mannucci 2019). Scaling relations, such as the mass-metallicity relation (MZR), provide key insights into the physical mechanisms involved in the growth and evolution of galaxies. At z ∼ 2 − 4, the MZR has been studied using both the gas-phase metallicity (e.g., Erb et al. 2006; Troncoso et al. 2014; Sanders et al. 2021) and the stellar metallicity (e.g., Sommariva et al. 2012; Cullen et al. 2019; Calabrò et al. 2021), resulting in the determination of a strong redshift evolution towards lower metallicities at a given stellar mass.

Moreover, as different chemical elements are produced by stellar populations on different timescales, the relative abundance of elements enables us to obtain constraints on the star formation history (SFH) of galaxies (Maiolino & Mannucci 2019). The C/O abundance ratio is a powerful indicator because most of the oxygen is synthesized in massive stars (> 10 M), while carbon is produced in massive and intermediate-mass stars. Thus, a time delay in the production of carbon and its ejection to the interstellar medium (ISM) makes C/O a measurable “chemical clock” for the relative ages of the stellar populations in galaxies (Garnett et al. 1995) and an important indicator for constraining chemical evolution models (Vincenzo & Kobayashi 2018).

The C/O ratio has been studied in local dwarf galaxies and HII regions of disk galaxies (e.g., Garnett et al. 1995; Chiappini et al. 2003; Esteban et al. 2014; Peña-Guerrero et al. 2017; Berg et al. 2019a). A continuous increase of C/O with O/H is found above one fifth solar metallicity, but the relation flattens at lower metallicities (12 + log(O/H) < 8), showing a significant scatter of C/O values for a given metallicity (Garnett et al. 1995; Berg et al. 2016; Pérez-Montero & Amorín 2017). A detailed comparison with models led Berg et al. (2019a) to conclude that the C/O ratio is very sensitive to the assumed SFH, in such a way that longer and lower star formation efficiency bursts lead to low C/O ratios. Chemical evolution models with different prescriptions have been developed to understand the evolution of C/O with metallicity (e.g., Carigi 1994; Henry et al. 2000; Chiappini et al. 2003; Mattsson 2010; Carigi & Peimbert 2011; Mollá et al. 2015; Vincenzo & Kobayashi 2018) but several variables remain unconstrained, especially at high redshifts. Therefore, measurements of C/O for galaxies at different metallicities remain crucial in the study of the properties of CIII] emitters and SFGs in general, given that they can lead to variations in their observed CIII] emission (e.g., Jaskot & Ravindranath 2016; Nakajima et al. 2018a).

In this work, we focus on the average properties of CIII] emitters using the spectral stacking technique. One key goal is to study, for the first time at this redshift, the relation between the mean stellar metallicity and C/O abundances of galaxies, which is discussed in terms of other physical properties of the sample. This will be useful for interpreting future observations with the James Webb Space Telescope (JWST) at higher redshifts, where only the rest-frame UV spectral lines would be accessible. In a forthcoming paper (Llerena et al., in prep.) we will present a second study based on individual galaxies.

This paper is organized as follows. In Sect. 2, we present the sample selection, the basic properties of the sample, and our stacking method. In Sect. 3, we present a qualitative and quantitative description of the emission and absorption lines detections via different emission-line diagnostics, the estimation of metallicities, C/O abundances, and different correlations found for our sample. In Sect. 4, we discuss our results, focusing on the stellar mass–metallicity relation and the C/O–metallicity relation. Finally, in Sect. 5, we present our conclusions.

Throughout the paper, we assume the following cosmology: ΩM = 0.3, ΩΛ = 0.7, H0 = 70 km s−1 Mpc−1. We adopt a Chabrier (2003) initial mass function (IMF). We consider the solar metallicity Z = 0.0142, log(O/H) = 8.69, and log(C/O) = −0.26 (Asplund et al. 2009). We assume, by convention, a positive EW for emission lines and rest-frame EW. We use the following notation for metallicity for consistency with local conventions: [X/Y] = log(X/Y)–log(X/Y).

2. Methodology

2.1. Sample selection

In this work, we use spectra from VANDELS (McLure et al. 2018; Pentericci et al. 2018), an ESO public spectroscopic survey conducted with VIMOS at the Very Large Telescope. VANDELS obtained unprecedented high S/N optical spectra of ∼2100 galaxies at redshift 1.0 ≤ z ≤ 7.0 in the UKIDSS Ultra Deep Survey (UDS: 02:17:38, −05:11:55) and the Chandra Deep Field South (CDFS: 03:32:30, −27:48:28) fields. Ultra deep spectra for every single galaxy has a minimum (maximum) total exposure time of 20 h (80 h), along with a mean spectral resolution R ∼ 580 and dispersion of 2.5 Å/pixel in the wavelength range of ∼4800−10 000 Å. The survey strategy and design, including target selection and data reduction, is fully described in a series of papers (McLure et al. 2018; Pentericci et al. 2018; Garilli et al. 2021). In short, VANDELS targets can be classified according to their selection criteria as bright SFGs in the range 2.4 ≤ z ≤ 5.5 and Lyman-break galaxies (LBGs) in the range of 3.0 ≤ z ≤ 7.0, as well as a smaller sample of passive galaxies (1.0 ≤ z ≤ 2.5) and AGN candidates. In this work, we only selected galaxies from the SFGs and LBGs targets.

Our sample is drawn from VANDELS DR3, which consists of 1774 galaxies – a subset of the 2087 galaxies included in the VANDELS final data release (Garilli et al. 2021). We selected galaxies with spectroscopic redshift quality flag 3 or 4, which means 95% and 100% of confidence in their spectroscopic redshift (McLure et al. 2018). We selected galaxies at 2 < z < 4 to ensure that the CIII] emission lines are included in the spectral range provided by the VANDELS spectra. The detection of CIII], typically the strongest nebular emission line in the sample, at S/N > 3 is required to ensure a proper measurement of the systemic redshift. With this constraint, from a parent sample of 746 galaxies with the above redshift range and quality flags, a first sample of 225 galaxies were selected based on their CIII] emission (130 in the CDFS field and 95 in the UDS field).

We cross-matched the sample of CIII] emitters with the 7 Ms CDF-S catalogue (Luo et al. 2017) and the ∼200−600 Ks X-UDS catalogue (Kocevski et al. 2018) in order to discard galaxies with X-ray emission within 3 arcsec of separation. We also discarded galaxies with spectral features consistent with AGNs or with strong sky residuals. A total of eight galaxies were excluded from the sample. A more detailed analysis on the AGN sample in VANDELS will be presented in Bongiorno et al. (in prep.). Our final sample of CIII] emitting galaxies is made of 217 galaxies (hereafter, the C3 sample), which represents ∼30% of the parent sample. Figure 1 shows the rest-frame spectrum of one of the galaxies in the C3 sample. Some of the expected UV absorption and emission lines are marked by vertical lines. Detected emission lines that are relevant to our study are shown in both in 1D and 2D spectra and marked in different zoom-in panels.

thumbnail Fig. 1.

Spectrum of UDS20394, one of the more intense CIII] emitters in the C3 sample, whose estimated parameters are log(M/M) = 9.32, SFR = 4.01 M yr−1, MFUV = −20, MKs = −20.46, EW(Lyα) = 19.2 Å, and EW(CIII]) = 12.3 Å. The green faint line is the de-redshifted VANDELS spectrum and the black line is the same but resampled by a factor of 2. The blue line in the upper panels is the error spectrum. The red line in the intermediate panels is the scaled sky spectrum.

2.2. Systemic redshift and basic properties of the sample

In order to prepare the C3 sample for the stacking procedure, we followed the methodology described in Marchi et al. (2019) to derive accurate systemic redshifts (zsys) using the nebular CIII] line. In Fig. 2, we present the resulting zsys distribution for the C3 sample, which spans the range of 2.17 − 3.82 (⟨zsys⟩ = 2.98, σ = 0.43). Compared with the spectroscopic redshifts (zspec) of the sample reported in McLure et al. (2018), the systemic redshifts are slightly larger with a mean difference Δ(zsys − zspec) = 0.002, which corresponds to ∼4 Å at the rest-wavelength of CIII].

thumbnail Fig. 2.

Systemic redshift distribution of the CIII] emitting galaxies in the C3 sample.

The physical properties of the C3 sample and the remaining galaxies in the VANDELS DR3 parent sample at 2 < z < 4 are obtained from the Spectral Energy Distribution (SED) fitting using the Bayesian Analysis of Galaxies for Physical Inference and Parameter EStimation (BAGPIPES1) code. BAGPIPES is a state-of-the-art Python code for modeling galaxy spectra and fitting spectroscopic and photometric observations (Carnall et al. 2018), which has now been applied to the VANDELS final data release (Garilli et al. 2021). For this paper, the BAGPIPES code is run by fixing the redshift and using the 2016 updated version of the Bruzual & Charlot (2003) models based on the MILES stellar spectral library (Falcón-Barroso et al. 2011) and updated stellar evolutionary tracks of Bressan et al. (2012) and Marigo et al. (2013). The stellar metallicity is fixed to 0.2 Solar and the nebular component is included in the model, assuming an ionization parameter log(U) = − 3. We chose to fix these parameters to typical mean values found in SFGs at similar redshift (e.g., Cullen et al. 2019; Runco et al. 2021) to minimize the effects of possible degeneracies affecting the models (see, e.g., Castellano et al. 2014). We note, however, that the results presented in subsequent sections remain unchanged even when we allow these parameters to vary within the typically observed ranges. Dust attenuation is modeled using the Salim et al. (2018) model. The SFH is parameterized using an exponentially increasing τ-model. We obtained a mean value for the timescale τ = 5.49 Gyr for the C3 sample (τ = 5.36 Gyr for the parent sample), which essentially implies constant star-formation. On the other hand, we obtained a mean age, namely, the time since the SFH began, of 228 Myr for the C3 sample (270 Myr for the parent sample). The ages obtained from the SED fitting are thus longer than the timescales (< 20 Myr) in which the EW(CIII]) changes with age, according to photoionization models assuming continuous star formation (Jaskot & Ravindranath 2016).

In Fig. 3, we present some relations between the parameters extracted from the SED fitting, such as stellar mass and rest-frame luminosity in different filters. The C3 sample is color-coded by SFR, which range from logSFR [M yr−1] = 0.13−2.89 (⟨logSFR [M yr−1]⟩ = 1.33, σ = 0.42). The stellar masses of the C3 sample span from log(M/M) = 8.54 to 10.40 with a mean value of ⟨log(M/M)⟩ = 9.41 (σ = 0.33). The FUV(1500) luminosity, tracing the young stellar component of galaxies, span between MFUV = −21.93 to −18.43 mag, with a mean value of ⟨MFUV⟩= − 20.55 mag (σ = 0.65). The Ks band, which better traces the evolved stellar component, ranges between −24.16 and −18.90 mag, with a mean value of ⟨MKs⟩= − 21.44 mag (σ = 0.85). As expected, the rest Ks-band luminosity is a good tracer of the stellar mass of the galaxies, showing little scatter in Fig. 3.

thumbnail Fig. 3.

Relations between the resulting BAGPIPES parameters in the VANDELS main-sequence galaxies. The galaxies of the C3 sample are color-coded by star formation rate. Gray points are VANDELS galaxies of the parent sample that were not included in the C3 sample.

Figure 4 demonstrates that the C3 sample is mainly located along the M–SFR main-sequence followed by the parent sample. Only a few galaxies at the higher stellar mass end (≳1010 M) appear offset to higher SFR. The C3 sample is therefore fairly representative of the VANDELS DR3 parent sample in this redshift range. The parameters shown in Fig. 3 are thus used for the stacking analysis of the global properties of the sample and their distributions are displayed in the histograms in Fig. 5.

thumbnail Fig. 4.

M–SFR relation with our sample color-coded by EW (CIII]). The orange solid line is the main sequence at z = 3, according to Speagle et al. (2014). Gray points are as in Fig. 3.

thumbnail Fig. 5.

Distribution of the resulting BAGPIPES parameters in the C3 sample. Left to right: stellar mass, FUV, and Ks band luminosity. The vertical dashed lines represent the ranges of each bin for stacking according to Table 1.

Besides the physical parameters obtained from BAGPIPES, we also measure the EW(CIII]) and EW(Lyα) in all galaxies of the C3 sample. For CIII], we use slinefit2 following a similar scheme to the one that will be explained in detail in Sect. 2.3 for the stacked spectra. The EW(CIII]) distribution is presented in the left panel in Fig. 6. The EW(CIII]) has a mean value of ⟨EW(CIII])⟩ = 3.98 Å (σ = 3.12 Å). While most galaxies show low EW(CIII]) < 5 Å, we find a small number of strong CIII] emitters with EW(CIII]) up to ∼20 Å (∼11% of the C3 sample with EW(CIII]) > 8 Å). The EW(CIII]) values are shown in the color-code of Fig. 4, where the M–SFR plane is shown. It can be noticed that the intense and faint CIII] emitters are above and below the main sequence, with some trends suggesting that the more intense CIII] emitter have lower stellar masses and then lower star formation rates.

thumbnail Fig. 6.

EW(CIII]) (left) and EW(Lyα) (right) distributions of the galaxies in the C3 sample. The vertical dashed lines are the limits of the bins used for stacking according to Table 2.

About half of our C3 sample is at z > 2.9 with Lyα observable in the spectral range. For these 105 galaxies, we use the EW(Lyα) obtained by Cullen et al. (2020). The distribution of such values is presented in Fig. 6. We find that the EW(Lyα) span a wide range from −48.37 to 99.79 Å, with a mean value of 12.69 Å (σ = 28.78). These values are significantly higher than the mean EW(Lyα)∼2 Å of the parent sample. Thus, the C3 sample includes both strong Lyα emitting galaxies and galaxies with weak or absent Lyα emission. About 34% of the C3 sample with Lyα included in the spectral range are considered Lyα emitters galaxies (LAEs, i.e., EW(Lyα) > 20 Å) consistent with what is commonly found for LGBs at these redshifts (e.g., Cassata et al. 2015; Ouchi et al. 2020). A smaller fraction of LAEs (∼23%) is found for galaxies with non-detections (S/N < 3) of CIII] in the parent sample.

2.3. Stacking procedure

In this paper, we are interested in the characterization of the mean physical properties for the CIII] emitters in VANDELS. For this reason, we performed a stacking analysis, binning the C3 sample by stellar mass and rest-frame luminosity and EW. This allows us to increase the S/N of the data and probe properties, such as stellar metallicity, which would not be possible in individual objects.

For this purpose, we separated the selected galaxies in five stellar mass bins of width ∼0.3 dex and six bins of luminosity of 0.5 dex each. This way, we have a significant number of galaxies per bin, as shown in Table 1 and marked by the vertical dashed lines in Fig. 5. As shown in Fig. 3, stellar mass and Ks luminosity are correlated, as the latter is a good tracer of the former. While we expect stacks in these two quantities to produce similar results, we decided to use both of them to evaluate any possible difference due to the larger dynamic range of the Ks luminosity.

Table 1.

Bin ranges used for the stack analysis with BAGPIPES parameters.

We also separated the sample in bins of EW(CIII]) and EW(Lyα). For the former, we separated them in three bins of 4 Å. And for Lyα, the C3 sample is restricted to the subsample of galaxies with zsys > 2.9 where Lyα is observable, either in absorption or emission. We chose bins of ∼20 Å for the EW(Lyα), which are presented in Table 2 and are marked by vertical dashed lines in Fig. 6.

Table 2.

Bin ranges used for the stack analysis with EW.

For the stacking, we use a non-weighted scheme following Marchi et al. (2017). All the individual spectra in the sample are first shifted into the rest-frame using the systemic redshift and then they are resampled onto a common grid according to the mean systemic redshift of the sample (zsys ∼ 2.98) and normalized to the mean flux between 1460 and 1540 Å. The final flux at each wavelength was taken as the median of all the individual flux values after a 3-σ clipping for rejecting outliers. The final wavelength range is where all spectra overlap and the spectral binning is 0.64 Å. The 1-σ error spectrum is estimated by a bootstrap re-sampling of the individual fluxes for each wavelength and the standard deviation of the resulting median stacked spectra. Changing the range of normalization to ∼1800 Å does not affect the shape of the stacked spectra.

We also tested alternative weighted schemes for stacking, similar to the one presented in Marchi et al. (2017) and used in Saxena et al. (2020), with a 1/σ2 weight, where σ is estimated as the flux error along the normalization range in each spectrum. The error spectra with the weighted scheme were larger compared with the median stacking. For this reason, we considered the median stacking for this work, as this offers a better representation of the global properties of the galaxies in each bin.

Four different median stacking schemes were performed depending on the redshifts included in each bin. In the remaining of this work, they are named as follows:

– Stack A: All the galaxies in the bin are stacked, considering the entire C3 sample. These stacks for each physical parameter are presented in Fig. 7.

thumbnail Fig. 7.

Resulting stacked spectra for each physical parameter with the C3 sample (what we call Stack A, see text for more details). Left to right: MKs, MFUV, and stellar mass. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line is the 1-σ error spectrum. The vertical lines mark known UV lines (in black: emission lines; in red: ISM absorption lines). Information about the number of galaxies, the mean redshift, and the mean parameter are included in each panel.

– Stack B: Only galaxies with z > 2.93 are considered for stacking in each bin. These galaxies have Lyα included in the spectral range. The stacks for each physical parameter are presented in Fig. B.1.

– Stack C: Only galaxies with z > 2.93 and with EW(Lyα) > 0 (i.e. Lyα in emission) are stacked. These stacks for each physical parameter are presented in Fig. B.2.

– Stack D: Only galaxies with z < 2.93 are stacked. In this subset, Lyα is not covered by the VANDELS spectra. Thus, we ignored whether these galaxies are Lyα emitters or not. The stacks for each physical parameter are displayed in Fig. B.3.

An additional subset of stacks A, B, C, and D by EW(CIII]) are presented in Fig. B.4. In the case of the stacks by EW(Lyα), only Stack B is performed. The resulting stacks are presented in Fig. 8.

thumbnail Fig. 8.

Resulting stacks by EW(Lyα). Details are the same as in Fig. 7. Information about the number of galaxies, the mean redshift, and the mean EW(Lyα) are included in each panel. Left panels: relative Lyα strength in each bin.

The above redshift dependence reduces the number of galaxies in each bin. In these cases, we adapt the appropriate binning for there to be at least four galaxies in each bin for stacking. This ensures the stack spectrum will gain at least a factor of 2 in the S/N. The final number of galaxies for each stacked spectrum is included in labels in Figs. 7, 8, and B.1B.4. We find the composite spectrum of a bin to be representative of the median properties of the galaxies in each bin. Small changes in the bin sizes used for the stacking may change the error bars in the derived parameters in the least populated bins, but they do not affect our results significantly. We discuss possible caveats related to the stacking analysis in Appendix C.

2.4. Line measurements

Emission-line fluxes and EWs of the lines in each stacked spectrum were measured using slinefit, a software capable of simultaneously measuring emission and absorption lines and the UV-NIR continuum. For this purpose, slinefit uses templates built with the Bruzual & Charlot (2003) stellar population models.

For our measurements, we included rest-frame UV emission and absorption lines at λrest > 1500 Å from Shapley et al. (2003). Rest-frame UV lines at λrest < 1500 Å were not included in the measurements. In particular, the Lyα line is instead measured following the same fitting technique presented in Cullen et al. (2020). In all the slinefit fitting runs, we allow for lines other than CIII] to have a small offset with respect to the systemic velocity and for the minimum width of the lines to be 100 km s−1.

For a first set of measurements, we only simultaneously measured CIII] and closer lines (AlIIIλ1855 Å, SiIIIλλ1883,1892 Å, MgIIλ2799 Å). We obtained the width of CIII] to be ∼300−350 km s−1 in all stacked spectra. We used this value to constrain the maximum width for the other emission lines.

Afterwards, for a second set of measurements, we measured HeIIλ1640 Å, OIII]λλ1666,1660 Å, and CIVλλ1548,1550 Å (hereafter, HeII, OIII], and CIV, respectively), with the constraint in the maximum width. In these cases, the maximum offset allowed is 100 km s−1, except for CIV for which a maximum of 1000 km s−1 is allowed because larger offsets are observed. The P-Cygni profile of CIV is fitted assuming the same intensity for both components. Both components of OIII] are fitted with their ratio unconstrained.

Due to the complex CIV profile, the measurements with slinefit are found to slightly underestimate the continuum. For this reason, a more detailed continuum determination was performed. First, the continuum is fitted with a linear function between 1400 and 2000 Å, masking out regions with emission and absorption lines detected. Then, the spectrum was continuum-subtracted. Finally, the CIV flux was found by direct integration of the emission line profile after imposing a maximum base line width of 4 Å (or ∼390 km s−1), that is, the typical value obtained for CIII]. The EW is estimated using the mean continuum flux in the same integrated range.

All the above measurements are performed for all the stacked spectra with the 0.6 Å/pixel sampling, but the S/N in the case of faint emission lines is low and then the measurements are additionally performed in the resampled spectra by a factor of 2. For the resampling, we use SpectRes3, a software that efficiently resamples the spectra and their associated uncertainties, preserving the integrated flux. An example of the measurements can be seen in Fig. 9. Henceforth, we considered the resampled spectra measurements for the emission lines, which are presented in Tables D.1D.5.

thumbnail Fig. 9.

Emission-line flux and EW measurements in the stack by EW(Lyα) > 20. Panel a: continuum-subtracted spectrum. The orange line shows the error spectrum (1σ). The blue shaded region are the pixel integrated for the CIV line. Panels b and c: red line is the slinefit fitting for the spectrum. The orange shaded region shows 1σ uncertainty of the stacked spectrum.

In the same tables, the mean color excess E(B − V) is reported for each stack. They are estimated from the individual E(B − V) values for each galaxy in each bin, which are obtained via a BAGPIPES fitting. In the C3 sample, E(B − V) range from ∼0.01−0.43 mag, with a mean value of 0.098 ± 0.014 mag. The mean E(B − V) of the stacks were used to compute the reddening correction4 using the Calzetti et al. (2000) extinction curve for simplicity and assuming that the color excess of the stellar continuum is the same as the color excess for the nebular gas emission lines. Despite the evidence that this assumption could not be true (e.g., Calzetti et al. 2000; Reddy et al. 2015) and that the ionized gas E(B − V) could be larger than a factor of 2.2 the stellar E(B − V) (in particular, galaxies with high SFR), we assumed this for simplicity. However, we note that the results of this paper are not affected if we change this prescription to more extreme assumptions. Using the calibration presented in Sanders et al. (2021) to correct the gas extinction from the SED extinction, we obtain a factor up to ∼3 of difference between gas and stellar extinction, but even with those values, the trends found in this paper are not altered significantly.

Line fluxes are presented uncorrected by extinction in Tables D.1D.5. However, the results and figures shown in subsequent sections are based on the dereddened quantities.

3. Results

3.1. Identification of UV absorption and emission lines

The high S/N spectra of the stacks allows us to identify several interesting features both in absorption and emission in the rest-frame UV spectra. Among these features, low-ionization interstellar lines such us SiIIλ1260 Å, OI+SiIIλ1303 Å, CIIλ1334 Å, SiIIλ1526 Å, FeIIλ1608 Å, and AlIIλ1670 Å are found. In Fig. 7, where the stacks A (redshift-independent) by stellar mass and luminosities are shown, we find that the stronger ISM absorption lines are in the stacks built with higher stellar masses (or more luminous at a given broadband). This is expected since these lines are saturated in low-resolution spectra and then their width increases with dynamical mass. The same trend is shown in Figs. B.1B.3 for the stacks B, C, and D by stellar mass and luminosities.

We find a similar trend when considering the stacks by EW(Lyα) in Fig. 8. We find the stronger low-ionization ISM absorption lines in the stacks with smaller EW(Lyα), namely, when Lyα is in absorption, while the ISM absorption lines are barely identified in the stack with the larger EW(Lyα). This is consistent with previous observations (e.g., Shapley et al. 2003). Regarding the stacks by EW(CIII]), we note the same trend as we show in Fig. B.4. The stronger low-ionization ISM absorption lines are found in the stacks with smaller EW(CIII]).

In addition to the low-ionization features associated with neutral outflowing gas, we identified high-ionization interstellar absorption lines such as SiIVλλ1393,1402 Å, CIV, and NVλλ1238,1242 Å. While SiIVλλ1393,1402 Å, and NVλλ1238,1242 Å are only identified in all the stacks B and C, and in the stacks by EW(Lyα) due to the spectral range, CIV is identified in all the stacks. We note in Figs. 7, 8, and B.1B.4 that the stronger absorption lines are in the stacks of higher stellar mass, brighter in luminosity, lower in EW(Lyα), and lower in EW(CIII]), which is in line with the trend observed in low-ionization ISM features.

In our stacks, we also identified fine-structure emission lines of SiII that have been observed in the rest-UV spectrum of star-forming galaxies (e.g., Shapley et al. 2003). SiII*λ1533 Å is in the spectral range of all the stacks. This faint line is between two ISM absorption lines and is identified in most of the stacks. This line is particularly more intense in the more massive and luminous galaxies showing both lower EW(Lyα) and lower EW(CIII]), but the trend is less clear than the one we find in the ISM absorption lines (see Figs. 7, 8, and B.1B.4). Other fine-structure lines, such as SiII*λ1265 Å and SiII*λ1309 Å, are also identified in stacks B and C, and in the stack by EW(Lyα) (see Figs. B.1, B.2, 8, and B.4). These lines are identified in all the stacks, irrespective of stellar mass, luminosity, or EW bin, suggesting they are a more common feature in the UV spectra of CIII] emitters.

In addition to CIII], we identified nebular emission lines such as OIII], CIV, and HeII, which are central to the main goals of this paper. We observe that the strength of the nebular lines depends on the stellar mass and luminosity. We find that in general the less massive (and fainter in any band) stacks show the more intense Lyα, CIV, HeII, OIII], and CIII] nebular lines. In the case of more massive (and brighter in any band), we find that the same set of nebular lines tend to be fainter. In particular for CIV and HeII, they show a stellar wind component as suggested by the P-Cygni profile (in the former) or broad profile (in the latter). In particular, CIV shows a P-Cygni type profile in all stacks with the emission being more intense in the less massive (or faintest) stacks and in the ones with lowest EWs, either CIII] or Lyα (see Figs. 7, 8, and B.1B.4).

In Fig. 8, it is worth noticing that the nebular features appear less dependent of the EW(Lyα) than the ISM absorption features, for which clearer differences are seen. We also note that upon comparing stacks B and C, which differ with regard to the inclusion (or not) of galaxies with Lyα in absorption, we find no strong difference in nebular emission lines or ISM absorption, but there are differences in the strength of Lyα.

3.2. Relation of EW(CIII]) with luminosity and stellar mass

In this section, we present a qualitative description of the relation of EW(CIII]) with the physical parameters used for stacking, which are shown in Fig. 10. We find a trend in which the stacks with more massive galaxies have lower EW(CIII]), while the more intense CIII] emitters correspond to the ones built with the lowest stellar mass; however, the scatter is large when we consider individual objects (small circles) and the relation is weak, especially for stacks A and D, which are those containing galaxies at z ≲ 3 (small green circles). Something similar is observed with the broad-band luminosities, where the stacks of fainter objects tend to have higher EW(CIII]). The scatter follows the one observed for the luminosity of the individual galaxies of the C3 sample. We note that galaxies in the C3 sample at z ≲ 3 (included in stack D) show a mean EW(CIII]) = 2.4 Å and tend to have lower EW(CIII]) than galaxies at z ≳ 3 (included in stack B and C), which show a mean EW(CIII]) = 5.3 Å. This is related to the overlap of VANDELS targets with slightly different selection criteria around z ∼ 3. For instance, galaxies selected as LBGs tend to show lower stellar masses than those selected as bright SFGs. We refer to Garilli et al. (2021) for a more detailed discussion on the effect of the VANDELS selection criteria on the galaxies physical properties.

thumbnail Fig. 10.

Relation of EW(CIII]) with (from left to right) stellar mass, Ks-band luminosity, and FUV luminosity, for the stacks A, B, C, and D for each parameter, and for the individual galaxies in the C3 sample (small circles; in green, we show galaxies at z ≲ 3, and in red, galaxies at z ≳ 3). In the x-axis, the parameters correspond to mean values for stacks, while for the C3 sample, the parameters correspond to the SED fitting values.

Comparing the different schemes of stacking, it can be noticed that the Stack B and C, that is, the stacks that only include galaxies at z ≳ 3, have larger EW(CIII]) than Stacks A and D for a given stellar mass or luminosity. On the other hand, Stack D, including only galaxies at z ≲ 3, tend to have lower EW(CIII]) for a given mass or luminosity than all the other stacking schemes. Except for the stack with the lowest stellar mass, EW(CIII]) ≳ 3 Å are not observed in Stack D. In addition, comparing Stack B and C, namely, the effect of including (or not) galaxies with Lyα in absorption, we note that slightly high EW(CIII]) are observed when only Lyα in emission is included.

A similar brief analysis is presented in Appendix A for the non-detected CIII] emitters from the parent sample. We find they are consistent with these results and they are intrinsically faint CIII] emitters with EW(CIII]) ≲ 2 Å for any FUV luminosity.

3.3. Diagnostic diagrams based on UV emission lines

Different diagnostic diagrams using rest-frame EWs of CIII], CIV, and OIII] and the line ratios of CIII], CIV, and HeII have been proposed to identify the main source of ionizing photons in distant galaxies (e.g., Nakajima et al. 2018a; Hirschmann et al. 2019). These diagnostics are useful for determining the nature of the dominant ionizing source of galaxies. However, they need to be constrained using large samples of galaxies where these lines are detected. Given that the C3 sample is draw from a parent sample of photometrically selected SFGs, they can be used to probe and constrain these models using nebular flux ratios and EWs of systems that are not particularly extreme in their properties.

In this subsection, we explore the ionization properties of the C3 sample using different diagnostic diagrams proposed by Nakajima et al. (2018a), which are presented in Figs. 11 and 12. These diagnostic diagrams consider equivalent widths and line ratios of UV nebular lines to classify star-forming and AGN-dominated galaxies, respectively. They are based on the prediction of photoionization models showing that UV lines are sensitive to the shape of the incident radiation field and can be used to distinguish the nature of the dominant ionizing source – that is, whether it comes from pure star formation or AGN (cyan and red circles in Fig. 12, respectively). We also include the results from synthetic emission lines from Hirschmann et al. (2019) that take into account composite galaxies classified by BPT diagram (Baldwin et al. 1981) for the limits in the diagnostic diagrams and the criteria for defining a composite galaxy depends on the ratio of black hole accretion rate and star-formation rate.

thumbnail Fig. 11.

Diagnostic diagrams with EW of UV emission lines for the different stack sets, color-coded by the parameter used for stacking. The black dashed lines separate AGNs and star-forming galaxies as proposed in Nakajima et al. (2018a). The dotted and dotted-dashed lines are from Hirschmann et al. (2019), meant to separate between SFGs and composite as well as composite and AGN, respectively. Top row: diagnostic of EW(CIII])–CIII]/HeII ratio. The gray shaded region is where the models overlap. Middle row: diagnosis of EW(CIV)–CIV/HeII ratio. Bottom row: diagnosis of EW(OIII])–OIII]/HeII ratio. The green boxes are the composites from Le Fèvre et al. (2019) and the higher the size, the higher the EW (CIII]). The red X mark is the composite from Amorín et al. (2017). The red rectangles are the stacks from Nakajima et al. (2018b) with the brightest MUV (smallest rectangle), smaller EW(Lyα), faintest MUV, and larger EW(Lyα) (largest rectangle).

thumbnail Fig. 12.

Diagnostic diagrams with fluxes of UV emission lines for the different stack sets, color-coded by the parameter used for stacking. Symbols are the same as in Fig. 11. The red and cyan points are the photoionization models for AGN (Feltre et al. 2016) and SF (Gutkin et al. 2016) for metallicity Z = 0.0001, 0.0002, 0.0005, 0.001, 0.002, and for ionization parameter log(U) = −4, −3, −2, −1. The dotted-dashed lines in the bottom panel separates composite and AGNs according to Hirschmann et al. (2019).

In Fig. 11, we present the following diagrams: EW(CIII]) versus CIII]/HeII, EW(CIV) versus CIV/HeII, and EW(OIII]) versus OIII]/HeII. Our results are color-coded by the physical parameter used for stacking and with symbols representing each type of stack (see Fig. 11). As a reference for comparison, we include in Figs. 11 and 12 similar results for stacks of CIII-emitters from the VUDS survey (Le Fèvre et al. 2015). The black dashed lines are the limits between star-forming galaxies (on the right of the lines) and AGNs (on the left of the lines). The gray shaded region in the first diagram is where the models overlap and the classification may be ambiguous. In stacks where one of the lines in the diagnostic ratios is not detected at 2σ, the 2σ-limits are taken into account. Instead, if two lines involved in a given diagnostic are undetected, then the stack is not considered.

Overall, the main result from Fig. 11 is that all the stacks explored in our VANDELS sample lay within the region dominated by ionization driven by star-formation. In the upper panel, stacks with lower stellar mass and fainter broad-band luminosity tend to show higher EW(CIII]) and higher CIII]/HeII ratios. Similar results are found in the middle panel of Fig. 11, which shows similar trend with EW(CIV) and CIV/HeII. Finally, the bottom panel of Fig. 11 shows the EW(OIII]) as a function of the OIII]/HeII ratio. Similar trends are found but with few stacks closer to the demarcation lines. Besides being consistent with the region dominated by star-formation, all the stacks lay in the region of composite galaxies according to Hirschmann et al. (2019) in all the diagnostic diagrams by EW.

We note that our stacks with stronger CIII] emission have consistent line ratios with regard to the stack of CIII] emitters of similar EW in VUDS (5 < EW(CIII]) < 10 Å, the smaller green rectangle in the top panel in Fig. 11) presented in Le Fèvre et al. (2019). The larger green symbols show VUDS galaxies with EW(CIII]) > 10 Å, which are closer to the demarcation lines. Therefore, our stacks explore a region in the emission line parameter space of photoionization models that is poorly explored in previous surveys and is occupied by the weak CIII] emitters.

In order to complement the above diagnostics with EW, we also explore diagnostic diagrams using only UV emission-line flux ratios. We compare them with the photoionization models from Feltre et al. (2016) and Gutkin et al. (2016) to probe the dominant ionizing source in our C3 sample. For the OIII] fluxes, we considered the sum of the OIII]λ1661 Å and OIII]λ1666 Å lines by adopting a theoretical ratio of OIII]1661 Å/OIII]1666 Å = 0.41 from Gutkin et al. (2016) for a log U = −2.

In the panels of Fig. 12 (from top to bottom), we present CIV/CIII] as a function of (CIV+CIII])/HeII ratio, and CIV/CIII], CIV/HeII and CIII]/HeII as a function of OIII]/HeII, respectively. Overall, the C3 stacks are fully consistent with diagnostics shown in Fig. 11, and point to pure star formation as the dominant source of ionization. Moreover, given that the galaxies in the C3 sample are selected to be star-forming and X-ray sources were excluded, these results lead to constrain photoionization models at z ∼ 3.

A few stacks in these diagnostic diagrams lie close to the demarcation lines. As shown in previous studies (e.g., Feltre et al. 2016; Gutkin et al. 2016; Nakajima et al. 2018a), these diagnostics may have some overlapping region in which both star-forming and AGN driven ionization may coexist. Only a few points lie at low OIII]/HeII and low CIII]/HeII (or CIV/HeII), a region where some contribution from AGN might be expected according to models. Also, Fig. 12 shows that the less massive galaxies tend to have higher OIII]/HeII ratios. In the bottom panel of Fig. 12, demarcation lines from Hirschmann et al. (2019) are shown. The stacks lay consistently in the region of composite galaxies.

As shown in Hirschmann et al. (2019), the demarcation lines between composite and SFGs are clearer at z = 0 − 1. At higher redshifts the models overlap in this zone. This could be an effect resulting from the evolution of the criteria for defining a composite galaxy with redshift or the change in the demarcation lines in the BPT diagram depending on the ionization parameter, electron density, or extreme UV ionization field (Kewley et al. 2019), which are likely to change due to the evolution of mass-metallicity relation with redshift. The higher electron temperatures at low metallicity can enhance the strengths of collisionally excited lines. It is also possible to find some pure SFGs with higher C/O in the composite region because the demarcation lines are based on fixed C/O at a given metallicity.

3.4. Apparent EW(Lyα)–EW(CIII]) correlation

Previous studies have reported a positive correlation between EW(Lyα) and EW(CIII]) (e.g., Shapley et al. 2003; Stark et al. 2014; Le Fèvre et al. 2019; Cullen et al. 2020). Such a relation is potentially useful for using CIII] to identify galaxies in the epoch of reionization where Lyα is strongly attenuated by the IGM.

In Fig. 13, we explore the possible correlation between EW(Lyα) and EW(CIII]) for our C3 sample. We consider spectra from Stack B by stellar mass, FUV, Ks luminosity, EW(CIII]), and EW(Lyα). We also include data from literature, both stacks (Shapley et al. 2003; Amorín et al. 2017; Nakajima et al. 2018b; Cullen et al. 2020; Feltre et al. 2020) and individual galaxies at similar redshifts (Stark et al. 2014), for which the existence of such correlation has been confirmed. We fit a linear regression to all the above data, which gives:

EW ( Ly α ) = ( 10.67 ± 0.89 ) × EW ( CIII ] ) ( 26.38 ± 5.72 ) . $$ \begin{aligned} \mathrm{EW}(\mathrm{Ly}\alpha ) = (10.67 \pm 0.89) \times \mathrm{EW(CIII]}) - (26.38 \pm 5.72). \end{aligned} $$(1)

thumbnail Fig. 13.

EW(Lyα)–EW(CIII]) relation. The results from stacks B are represented by the triangle symbols, with different color depending on the parameter used for stacking: blue for stellar mass, yellow for FUV luminosity, magenta for Ks luminosity, red for EW(CIII]), and green for EW(Lyα). Individual galaxies with measured Lyα from the C3 sample are the small blue circles with their typical errors on the upper left. Previous results from literature at similar redshift are also displayed from stacking (Shapley et al. 2003; Amorín et al. 2017; Nakajima et al. 2018b; Cullen et al. 2020; Feltre et al. 2020) and individual objects (Stark et al. 2014). Dashed black line is the best fit in Eq. (1) including our stacks and the sample from literature. Dashed red line is the best fit in Eq. (2), including just our stacks.

The linear fit is performed using LMFIT (Newville et al. 2014) with a least-squared method weighted by 1/σ, where σ is the uncertainty in the parameter used for fitting. We use the same method for linear fitting in the following sections. When the linear fit is displayed in the figures, the shaded region is the 3-σ uncertainty band of the fitting and it covers the range where the fit is performed and valid.

The relation in Eq. (1) fits the data with a relatively large scatter, which is shown at 3σ level by the grey band in Fig. 13. If we only include our stacks in the linear fit, we find that the best fit is

EW ( Ly α ) = ( 2.92 ± 0.85 ) × EW ( CIII ] ) ( 4.65 ± 2.36 ) , $$ \begin{aligned} \mathrm{EW}(\mathrm{Ly}\alpha ) = (2.92 \pm 0.85) \times \mathrm{EW(CIII]}) - (4.65 \pm 2.36), \end{aligned} $$(2)

which has a lower slope. The shallower relation can be an effect for the low EW(CIII]) of our stacks and because most of the literature sample is selected by Lyα or by strong CIII] emission. If we compare our stacks built by EW(Lyα) (green symbols in Fig. 13), we find that they are in better agreement with Eq. (1).

Our VANDELS data allow us to probe the low EW end of this relation. Figure 13 shows that there is a clear trend which seems to hold even in the absence of very strong CIII] emitters in our stacks. We caution, however, that the functional form in Eqs. (1) and (2) is representative for the average population of star-forming galaxies at z ∼ 2 − 4 and that strong deviations for individual galaxies can exist, especially at low EW. This is illustrated in Fig. 13, where we show the distribution of individual galaxies in the C3 sample. We note that these objects are not included in the linear fit in Eqs. (1) and (2). The scatter shown by the individual galaxies is larger than the typical uncertainties. A similar diagram was presented in Marchi et al. (2019) for a number of individual VANDELS sources.

Based on Cloudy photoionization models, Jaskot & Ravindranath (2016) show that at a given EW(Lyα), the scatter in EW(CIII]) can be as high as 10−20 Å, which is comparable to the observed level of scatter among galaxies, depending on the different metallicities, ionization parameters, and ages considered for the models. In general, higher EW(CIII]) for a given EW(Lyα) indicates a higher ionization parameter, younger ages, and lower stellar metallicity. In a future work, we will address this point by using our large sample of CIII] emitters to constrain different photoionization models with their individual measurements.

Overall, the correlation found in Fig. 13 suggests that CIII] emitters are good markers of LAEs, especially for galaxies with low stellar mass, low luminosity, and high star formation rates. This confirms the potential use of CIII] to identify and study galaxies at the epoch of reionization, for which Lyα emission is strongly attenuated due to IGM opacity. However, this could be challenging due to the lower EWs of CIII] compared to that of Lyα in SFGs. A similar stacking approach will be useful with large enough samples at z > 6 in studies of their global properties, but CIII] may be the only robust and high S/N emission that may be observed to have individual detections.

3.5. Relation between stellar metallicity and the CIII] and Lyα equivalent widths

The rest-frame UV spectrum is dominated by the continuum light from young, massive stars that contains features of the chemistry of stellar photospheres and expanding stellar winds. The strength of these features have a strong dependence on the total photospheric metallicity. The full UV-spectrum fitting uses the strength of these faint photospheric features to estimate the stellar metallicity.

We study the stellar metallicity (tracing the Fe/H abundance) of the C3 sample following two complementary approaches. First, we follow the full spectrum fitting method described in Cullen et al. (2019). In short, the high S/N stacked spectra are fitted using a Bayesian approach with Starburst99 (SB99) high-resolution WM-Basic theoretical models with constant star-formation rates (Leitherer et al. 2014). As a result, for some stacked spectra the stellar metallicity is an upper limit because the model parameter space does not extend below Z = 0.001 (0.07 Z). In these cases, a 2-σ upper limit is reported.

A second set of stellar metallicity estimates is performed using the method presented in Calabrò et al. (2021). This method is based on stellar photospheric absorption features at 1501 Å and 1719 Å, which are calibrated with SB99 models and are largely unaffected by stellar age, dust, IMF, nebular continuum, or interstellar absorption. These estimations were only possible in stacks A and B, which are the ones with the highest S/N (∼20−30). Using this method, we find consistent stellar metallicity values with those obtained with the full spectral fitting. However, the stellar metallicities based on the two photospheric indices show larger uncertainties (up to ∼0.6 dex). Hereafter, we use the results from the first approach.

The stellar metallicity of the C3 sample ranges from log(Z/Z) = − 1.09 (∼8% solar) to −0.38 (∼40% solar), with a mean value of log(Z/Z) = − 0.8 (∼16% solar). All stellar metallicities for stacks are reported in Tables D.1D.5.

We explore the relation of the stellar metallicity with the EW(CIII]), which is shown in Fig. 14. For this, we use the spectra from Stack A, which include the entire C3 sample of galaxies, irrespective of redshift. We include the stacks by stellar mass, FUV luminosity, Ks luminosity, and EW(CIII]) for the fitting. The best linear fit to data is:

log ( Z / Z ) = ( 0.51 ± 0.08 ) × log ( EW ( CIII ] ) ) ( 0.57 ± 0.03 ) , $$ \begin{aligned} \log (Z_{\star }/Z_{\odot }) = (-0.51 \pm 0.08) \times \log (\mathrm{EW(CIII])}) - (0.57 \pm 0.03), \end{aligned} $$(3)

thumbnail Fig. 14.

EW(CIII])-stellar metallicity relation. The results from our A stacks are represented by star symbols. Colors for our stacks are the same as in Fig. 13. Black dashed line is the best fit in Eq. (3) and the gray shaded region is the 3-σ band uncertainty.

with EW(CIII]) in Å. We find a decrease of EW(CIII]) with stellar metallicity. Comparing these results with Cullen et al. (2020), we find an offset towards higher EW(CIII]) for a given stellar metallicity. We find three reasons for such offset. Firstly, the two samples only overlap in a narrow redshift range, as the Cullen et al. (2020) sample include galaxies between z = 3 and z = 5, thus excluding galaxies at z < 3 which are numerous in our sample. Secondly, the selection criteria for our sample is based on CIII], leading to stacks with higher EW(CIII]). Finally, the use of accurate systemic redshifts based on CIII] line profile fitting lead to stacks with higher EW(CIII]) than stacks for the entire parent sample using the VANDELS spectroscopic redshift. We find the latter can explain up to a difference of log(EW(CIII])) ∼ 0.2 dex in the stacks of lower stellar mass and luminosity.

In Fig. 15, we show stellar metallicities as a function of EW(Lyα) for the set of spectra from Stack B using different colors for stacks by stellar mass, EW(CIII]), and FUV, Ks luminosity, and EW(Lyα). We perform a linear fitting to only our stack data, excluding two bins for which the EW(Lyα) is negative (i.e., Lyα is in absorption). Our best fit is

log ( Z / Z ) = ( 0.30 ± 0.07 ) × log ( EW ( Ly α ) ) ( 0.58 ± 0.06 ) , $$ \begin{aligned} \log (Z_{\star }/Z_{\odot }) = (-0.30 \pm 0.07) \times \log (\mathrm{EW}(\mathrm{Ly}\alpha )) - (0.58 \pm 0.06), \end{aligned} $$(4)

thumbnail Fig. 15.

EW(Lyα)-stellar metallicity relation. The symbols and colors for our stacks are the same as in Fig. 13. Black dashed line is the best fit to the stacked data (triangles) presented in Eq. (4), with the gray shaded regions at 3-σ, respectively. The black symbols are the stacks in Cullen et al. (2020) with Lyα in emission.

with EW(Lyα) in Å. We find a decrease of EW(Lyα) with stellar metallicity. The relation in Eq. (4) provide larger metallicities at fixed EW(Lyα) compared to the one presented by Cullen et al. (2020) which is based on stacking of galaxies at z > 3 and include spectra with Lyα in absorption. In our fit, instead, we only consider stacks with Lyα in emission. We note, however, that despite the difference in redshift, the two stacks from Cullen et al. (2020) (black rectangles in Fig. 15 and built by binning in EW(Lyα)) with Lyα in emission show a trend that is are consistent with our stacks based on EW(Lyα) (green triangles).

Overall, the relation found in Fig. 15 confirms and adds robustness to the anticorrelation found for VANDELS galaxies out to z ∼ 5 in Cullen et al. (2020). We demonstrate in Figs. 14 and 15 that galaxies with stronger CIII] emission show larger Lyα EW and lower stellar metallicities of ≲10% solar. In correlations involving Lyα, it is important to note that the Lyα emission is resonantly scattered and the correlations are not easily interpreted as they typically are for the CIII] or other nebular emission lines. As shown in Cullen et al. (2020), Eq. (4) indicates that harder ionizing continuum spectra emitted by low metallicity stellar populations plays a role in modulating the Lyα emission in star-forming galaxies.

3.6. C/O ratio and its relation with EW and physical parameters

The relative abundances of carbon, nitrogen and other alpha elements to oxygen may provide insight not only into the origin of carbon in galaxies but also in their chemical evolution. However, constraining the C/O abundance is often difficult. For local galaxies, the emission lines often used to derive C/O are exceedingly faint carbon recombination lines (e.g., Esteban et al. 2014) or the CIII] collisionally excited line, which is accessible for low-metallicity objects only from space (e.g., Senchyna et al. 2017; Berg et al. 2019a). At z ≳ 1, the required emission lines lie in the optical range but even for low-metallicity objects their faintness require very deep observations or stacking (e.g., Shapley et al. 2003; Amorín et al. 2017). While this makes the C/O difficult to constrain, this abundance ratio is essential to understand different emission-line diagnostics (e.g., Feltre et al. 2016; Jaskot & Ravindranath 2016; Nakajima et al. 2018a; Byler et al. 2018) and, more generally, the origin of carbon and the chemical evolution of star-forming galaxies. Here, we explore the C/O ratio, which can be derived from the observed C and O lines in the UV.

We estimate the C/O abundance using the code HII-CHI-mistry in its version for the UV (Pérez-Montero & Amorín 2017) (hereafter, HCM-UV5) considering PopStar stellar atmospheres (Mollá et al. 2009) for the photoionization models used by the code. This PYTHON code derives the carbon-to-oxygen ratio (i.e., log(C/O)) from a set of observed UV emission-line intensities, which are also used to estimate ionization parameter and gas metallicity in a consistent framework with results provided by the direct Te-method. More details on this methodology can be found in Pérez-Montero & Amorín (2017). For C/O, we use as an input the CIV, CIII], OIII] fluxes (and their errors), which are reported in Tables D.1D.3, after extinction correction (as explained in Sect. 2.4). In most of the stacks, the OIII]λ1660 Å is not detected at 3σ. For this reason, we consider the theoretical ratio of OIII]λ1660 Å/OIII]λ1666 Å ∼ 0.4 from photoionization models with an ionization parameter of −3 and −2 (Gutkin et al. 2016). In the code, 25 Monte Carlo iterations are performed to estimate the uncertainties in the C/O calculations. In the cases where OIII]λ1666 Å is detected with a < 2σ level or when the C/O uncertainties are larger than 0.9 dex, C/O is estimated as a lower limit. Results are reported in Tables D.1D.5. We find log(C/O) values ranging from −0.68 (38% solar) to −0.06 (∼150% solar) with a mean value of −0.50 (60% solar).

Alternatively to HCM-UV, we used the empirical calibration between C3O3 (≡log((CIII] + CIV)/OIII])) and C/O found by Pérez-Montero & Amorín (2017) based on a control sample with C/O and metallicities obtained from UV and optical lines. Using this calibration log(C/O) = −1.07 + 0.80 × C3O3, which essentially provides an accurate fit to models predictions for the C3O3 index, we find consistent results, within the typical error of ∼0.2 dex, with those of HCM-UV. Small differences can be attributed to slight changes in the ionization parameter (see Fig. 2 in Pérez-Montero & Amorín 2017), which is constrained by HCM-UV using CIV and CIII]. All C/O results presented in subsequent analysis and figures are derived with HCM-UV but they are fully consistent with the C3O3 calibration.

We explore the relation of log(C/O) with EW(CIII]) in Fig. 16. We only include the results from stack A by stellar mass, luminosities and EW(CIII]), and the stacks by EW(Lyα). A linear regression to data gives

log ( C / O ) = ( 0.19 ± 0.14 ) × log ( EW ( CIII ] ) ) ( 0.36 ± 0.07 ) , $$ \begin{aligned} \log \mathrm{(C/O)} = -(0.19 \pm 0.14) \times \log (\mathrm{EW(CIII])}) - (0.36 \pm 0.07), \end{aligned} $$(5)

thumbnail Fig. 16.

Relation between the EW(CIII]) and C/O ratio as derived using HCM-UV. The red dashed line corresponds to the solar value. The cyan shaded region is the average C/O for Lyman-break galaxies in Shapley et al. (2003) and the red circle is the result for a composite spectrum of CIII]-emitters in Amorín et al. (2017). Individual CIII] emitters at high-z (red squares, Amorín et al. 2017) and low-z (black squares and circles, Berg et al. 2019a; Senchyna et al. 2021) are also included. Colors for our stacks are the same as in Fig. 13. The black dashed line is the best linear fit to our stacks, as displayed in Eq. (5).

with a Pearson correlation coefficient of rp = −0.30. This relation suggests a decrease of C/O abundance ratio with EW(CIII]), namely, the more extreme CIII] emitters tend to have lower log(C/O). We show in Eq. (3) that strong CIII] emitters also have low stellar metallicities, which lead to less cooling and higher nebular temperatures that enhance the CIII] emission. Therefore, Eq. (5) suggests a change in stellar metallicity. The relation between C/O and stellar and gas metallicities is discussed in Sect. 4.2.

We also observe a weak relation (rp = −0.29) between C/O with EW(Lyα), where LAEs tend to have lower C/O, thus suggesting higher EW(CIII]) and lower metallicity, namely, a younger chemical age. For EW(Lyα) = 20 Å, a log(C/O) ∼ −0.5 is found (∼60% solar) and a corresponding EW(CIII]) ∼ 10 Å (from Eq. (5)). We stress, however, that the relation between both C/O and EW(Lyα) might not be physically motivated and it relies on previous correlations found between EWs in Eq. (2).

On the other hand, in Fig. 17 we present the results of C/O for the C3 sample as a function of the physical parameter used for stacking and color-coded by EW(CIII]). We find an apparent mild increase of C/O with stellar mass (left panel on Fig. 17). A linear fit to data after excluding lower limits gives:

log ( C / O ) = ( 0.03 ± 0.12 ) × log ( M / M ) ( 0.72 ± 1.20 ) , $$ \begin{aligned} \log \mathrm{(C/O)} = (0.03 \pm 0.12) \times \log (M_{\star }/M_{\odot }) - (0.72 \pm 1.20), \end{aligned} $$(6)

thumbnail Fig. 17.

Relation between the different physical properties used for the stacks and C/O ratio using as derived using HCM-UV. The red dashed line corresponds to the solar value. The cyan shaded region is from Shapley et al. (2003). The results of the stacks are color coded by EW(CIII]). The black dashed line is the best linear fit for each parameter in Eqs. (6)–(8), respectively.

with rp = 0.07. This relation is weak and of limited use due to the lack of reliable C/O estimations in the high mass end for which OIII] is barely detected in our stacks. This makes the dynamical range of stellar mass too small to provide a more robust relation.

In the Ks band (middle panel on Fig. 17), we perform a linear fitting excluding the lower limits, which gives:

log ( C / O ) = ( 0.03 ± 0.03 ) × M K s ( 1.21 ± 0.72 ) , $$ \begin{aligned} \log \mathrm{(C/O)} = -(0.03 \pm 0.03) \times M_{K_s} - (1.21 \pm 0.72), \end{aligned} $$(7)

with rp = −0.37. We find an increase of C/O with Ks luminosity, but again the relation is weak mostly due to the lower limits for high-luminosity stacks.

On the other hand, the FUV luminosity (right panel of Fig. 17) has a stronger correlation with C/O. A linear regression excluding lower limits gives:

log ( C / O ) = ( 0.14 ± 0.03 ) × M FUV ( 3.47 ± 0.63 ) , $$ \begin{aligned} \log \mathrm{(C/O)} = -(0.14 \pm 0.03) \times M_{\rm FUV} - (3.47 \pm 0.63), \end{aligned} $$(8)

with rp = −0.74. This correlation clearly shows an increase of C/O in galaxies with higher FUV luminosity and can be used to estimate a mean C/O value from a galaxy luminosity. Assuming the FUV luminosity is a tracer of the recent SFR and that more evolved stellar populations may have a larger contribution in the C/O relation with stellar mass and Ks luminosity, the above differences might be explained invoking a strong dependence of C/O with star formation histories. However, larger samples and deeper spectra, especially for high-mass galaxies, would be needed to provide better insight.

We note that our results suggest that stacking by FUV luminosity results in a more homogeneous distribution of C/O within the bins, compared to the stellar mass and Ks luminosity. In the short dynamical range for stellar mass (∼1 dex where C/O was estimated), there is no clear trend for this sample.

In summary, the C/O abundance of CIII] emitters increases from less than half solar for fainter (low-mass) galaxies to about solar abundance for our brighter (high-mass) objects. The stronger CIII] emitters, namely, higher EWs, are found in low- luminosity, low-C/O galaxies.

3.7. Relation of SFR with EW(CIII]) and C/O

Here, we explore the relation between EW(CIII]) and C/O abundance ratio with SFR. In both cases, we illustrate these findings using the sample of Stack A by stellar mass, FUV and Ks luminosities, and EW(CIII]). We estimated the SFR with the mean value of the SED-based SFR of individual galaxies in each bin. In Fig. 18, we present these results. Overall, we find that EW(CIII]) decreases with SFR, as we also find for individual galaxies in the C3 sample. This is consistent with the trend observed in Maseda et al. (2017) for a smaller sample of strong CIII] emitters. Our stacks show that galaxies with SFR ≳ 30 M yr−1 have EW(CIII]) ≲ 3 Å. Instead, for galaxies with SFR ≲ 30 M yr−1, the relation becomes steeper and shows more dispersion on the EW(CIII]) values, with the more intense CIII] emitters having lower SFRs. Stacks with average SFR ≲ 10 M yr−1 show all EW(CIII]) larger than ∼3 Å. We note that the result in Fig. 18 does not imply that higher SFR tends to suppress CIII] emission, as this relation hides an underlying metallicity dependence that is key for interpreting the CIII] emission (see Sect. 4.1 for a discussion on stellar mass–metallicity relation). Higher SFR implies a larger number of ionizing photons, which tend to enhance CIII] emission. However, since our sample lies in the star-forming main sequence, galaxies with higher SFR have higher stellar mass and also higher metallicity, which implies more efficient cooling and weaker collisionally-excited lines such as CIII].

thumbnail Fig. 18.

Relation between SFR–EW(CIII]) (top panel) and SFR–C/O (bottom panel). Stars symbols are Stack A color-coded depending on the parameter used for stacking, according to legend, as in Fig. 13. Top panel: the black squares are results from CIII] emitters in Maseda et al. (2017). The small dots are the C3 sample with same colors as in Fig. 10. Bottom panel: the dashed black line is the best fit for the SFR–CO relation in Eq. (9).

Finally, we do not find a correlation of EW(CIII]) with specific SFR (sSFR = SFR/M). For the stacks shown in Fig. 18 and based on the mean values of stellar mass and SFR from the SED fitting, we find values for log(sSFR [yr−1]) ranging between −8.3 and −7.9 with no clear correlation with EW(CIII]). We note, however, that the dynamical range in sSFR probed by the C3 sample appears too small compared with the uncertainties to probe possible correlations with other observables.

On the other hand, Fig. 18 shows a trend between SFR and C/O abundance. The C/O tends to increase with the average SFR. This relation shows a large scatter and appears mainly driven by the stacks by MFUV. This result shows consistency with the relation found between C/O and MFUV in Fig. 17, as MFUV is a good tracer of SFR. A linear fit to this relation, excluding upper limits, gives

log ( C / O ) = ( 0.12 ± 0.11 ) × log ( SFR SED ) ( 0.61 ± 0.13 ) , $$ \begin{aligned} \log \mathrm{(C/O)} = (0.12 \pm 0.11) \times \log (\mathrm{SFR_{SED}})-(0.61 \pm 0.13), \end{aligned} $$(9)

where SFR is in M yr−1 (rp = 0.51). Considering that our C3 sample is representative of main-sequence galaxies, this relation shows the limitation of Eq. (6) because faint OIII] lines are not detected in stacks with high-mass galaxies.

It is worth noting that in Figs. 1416, and 18, different dynamic ranges in stellar metallicity and C/O are obtained depending on the parameter used for stacking. Overall, the trends found in these and other relations are mainly driven by the stacks binned by luminosity, which are more uniformly populated. In Appendix C, we show that our results remain unchanged if we restrict the fitting to only those stacks binned by luminosity.

4. Discussion

4.1. Stellar mass-metallicity relation of CIII] emitters at z ∼ 3

Scaling relations, such as the stellar mass-metallicity relation (MZR), are important diagnostics for understanding the evolution of galaxies. In particular, the MZR is shaped by different physical processes such as strong outflows produced by stellar feedback, infall of metal-poor gas, the so-called stellar mass “downsizing” – for which high-mass galaxies evolve more rapidly and at higher redshifts than low-mass ones – or by the shape of the high-mass end of the IMF (see Maiolino & Mannucci 2019, for a review). Thus, the MZR of a galaxy population in a determined redshift may provide clues on the dominant processes that affect its evolution in that period.

The gas-phase MZR (hereafter MZRg) has been explored in the local universe, exploiting the methods based on optical emission lines applied to large samples of galaxies (e.g., Tremonti et al. 2004; Andrews & Martini 2013). In the same way, the MZRg has been explored at higher redshifts (e.g., Erb et al. 2006; Mannucci et al. 2010; Pérez-Montero et al. 2013; Troncoso et al. 2014; Lian et al. 2015; Sanders et al. 2021), thus providing a clear evolutionary picture up to z ∼ 3.5 (e.g., Maiolino et al. 2008; Mannucci et al. 2010). The MZRg is found to evolve with redshift, with metallicity declining with redshift at a given stellar mass. In one of the latest works, Sanders et al. (2021) use Te-consistent metallicity calibrations to derive O/H for star-forming galaxies out to z ∼ 3.3 and explore the MZRg evolution. They find similar slopes at all redshifts for log(M/M)≲10 and a nearly constant offset of about 0.2−0.3 dex towards lower metallicities, as compared to local galaxies at a given stellar mass. After comparing with chemical evolution models, the authors argue that this is driven by both higher gas fraction (leading to stronger dilution of ISM metals) and higher metal removal efficiency, for instance, via feedback.

On the other hand, studying the redshift evolution of the stellar MZR (hereafter MZRs) has historically been more challenging due to the required high S/N continuum spectra. In the local universe, first studies of the MZRs by Gallazzi et al. (2005) and subsequent work via stacking of SDSS optical spectra for statistical samples of galaxies (Zahid et al. 2017), found stellar metallicity increasing over a large range of log(M/M)∼9 − 11. At higher redshifts, however, the lack of high S/N optical spectra for statistical samples precludes similar analyses. Estimates of stellar metallicity are found, instead, from metallicity-sensitive indices or full spectral modeling of deep rest-frame UV spectra sampling young, massive stars (e.g., Sommariva et al. 2012; Cullen et al. 2019; Calabrò et al. 2021). Therefore, they provide Z values expected to be similar to those derived for the ISM out of which young stars have recently formed. Recent studies have shown an evolving MZRs up to z ∼ 3, decreasing Z at fixed M by more than a factor of 2 from z = 0 to z = 3.5, similarly to what is found for the MZRg (e.g., Cullen et al. 2019, 2020; Calabrò et al. 2021).

Here, we explore the position of normal galaxies with detected CIII] emission in the MZR at z ∼ 2 − 4, while comparing with some of the above results from the literature. Using the results presented in previous sections, we probe the stellar mass-metallicity relation MZRs on the left panel in Fig. 19. In this figure, we use stacks A, B, and D computed by stellar mass. Our best linear fit to these points gives the following relation:

log ( Z / Z ) = ( 0.33 ± 0.10 ) × [ log ( M / M ) 10 ] ( 0.54 ± 0.05 ) , $$ \begin{aligned} \log (Z_{\star }/Z_{\odot }) = (0.33 \pm 0.10) \times [\log (M_{\star }/M_{\odot })-10] - (0.54 \pm 0.05), \end{aligned} $$(10)

thumbnail Fig. 19.

Relations between stellar mass and metallicity. Left panel: MZRs with the stacks A, B, D by stellar mass with star, triangle and circle symbols, respectively, and color-coded by EW(CIII]). The black symbols are stacks from Cullen et al. (2019, 2021), Calabrò et al. (2021). The orange circles is the local MZRs from Zahid et al. (2017). The blue dashed line is the best fit with the shaded region at 3-σ in Eq. (10), and the black dashed line is the MZRs from Cullen et al. (2019). Right panel: MZRg with the same stacks as left panel but re-scaled with the α-enhancement. Green dashed line from Steidel et al. (2014). Red dashed line from Troncoso et al. (2014). Black symbols are values from literature (Sanders et al. 2021; Cullen et al. 2021). The orange points are the local MZRg from Andrews & Martini (2013). The Zg scale in the y-axes corresponds to the gas-phase metallicity obtained from Eq. (11) (see text for more details).

which is represented by the blue dashed line in Fig. 19. The MZRs of CIII] emitters at 2 < z < 4 have a nearly constant offset of ∼0.4 dex compared with the local MZRs from Zahid et al. (2017) over more than one decade in M. Compared to other MZRs at similar redshift from Cullen et al. (2019, 2020) and Calabrò et al. (2021), we find a relatively good agreement in slope and normalization – thus suggesting that CIII] emitters are not different from their parent sample of star-forming galaxies in the MZRs. While stacks A and D are slightly shifted to high Z, stacks B appear more consistent with the MZRs derived in the above previous VANDELS studies. These small differences may arise from the different selection criteria. In particular, redshift selection used in these works only include galaxies at z > 3, while our D stacks include only galaxies at z < 2.9 and stacks A include a mix of galaxies below and above z = 2.9. We also find that for a given stellar mass, the intense CIII] emitters tend to lie below the MZRs, while the faint emitters, tend to lie above the MZRs (see Fig. 19). On the other hand, the offset between Eq. (10) and the MZRs reported in Cullen et al. (2020) is explained in part due to the difference in the redshift range covered by the samples and differences in the assumptions leading to stellar mass derivations. In Cullen et al. (2020), stellar masses are derived from SED fitting assuming solar metallicity and models that do not account for nebular emission. These assumptions lead to an offset of around ∼0.2 dex towards higher stellar masses compared with our updated catalog.

If we assume that stellar and gas-phase metallicity are the same, our results are consistent with the MZRg derived from Troncoso et al. (2014) (at 1-σ) for the range of stellar mass covered by our stacks. However, this assumption is not necessarily true, especially at high-z.

Sommariva et al. (2012) found a small difference (∼0.16 dex) between Zg and Z, but their conclusion was not robust given their large reported uncertainties. One might expect small differences because we are measuring stellar metallicity using UV absorption features driven by massive stars with short lifetimes and similar properties of the interstellar gas where they were formed. However, larger differences can be found for galaxies whose ISM has been enriched primarily by core collapse supernovae with highly super-solar O/Fe, as discussed in Steidel et al. (2016), Topping et al. (2020). Such conclusion has been recently reached in a work by Cullen et al. (2021), where it was found that a subset of galaxies in VANDELS at z ∼ 3.4 are α-enhanced (i.e., their O/Fe ratios are more than two times solar) from a direct comparison of their stellar and gas metallicities.

Studies of the MZRg using exclusively the rest-UV spectrum has strong limitations due to the lack of hydrogen lines besides Lyα. While this line has been used in galaxies with extremely high EWs (Amorín et al. 2017), we avoid the use of Lyα in our VANDELS sample because it is generally affected by resonant scattering and absorption by the IGM. An alternative method to constraining gas-phase metallicity is using nebular HeII instead of Lyα. This possibility has been implemented, for instance, in HCM-UV (version 4) and in Byler et al. (2020) using a He2-O3C3 calibration. However, HeII is generally weak in most stacks and it may include both nebular and stellar origin, which being difficult to disentangle is thus an additional source of error.

We estimate gas-phase metallicities for our stacks using the above two methods. When we compare them with the derived stellar metallicities we find a mean difference of log Z − log Zg ∼ 0.16 dex for the stacks by stellar mass, FUV luminosity, and EW(CIII]). However, the dispersion is larger (Δ(log Z − log Zg)∼0.25 dex) and since the uncertainties for the metallicities are also larger (up to ∼0.4 dex), the above comparison does not provide a robust assessment of the true difference between stellar and gas-phase metallicities, which are proxies of the Fe/H and O/H abundances, respectively.

For this reason, we follow an alternative approach to estimate gas-phase metallicities in our stacks. Following Cullen et al. (2021), we consider that

log ( Z g / Z ) = log ( Z g / Z ) + log ( Z / Z ) , $$ \begin{aligned} \log (Z_g/Z_{\odot }) = \log (Z_g/Z_{\star }) + \log (Z_{\star }/Z_{\odot }), \end{aligned} $$(11)

where log(Zg/Z)∼[O/Fe] is a proxy of the α-enhancement, which depends on stellar mass. Then, adopting the difference found by Cullen et al. (2021) for MZRs and MZRg, we infer that [O/Fe] ∼ 0.37−0.40 dex for the range of stellar masses of our sample. Thus, we use [O/Fe] ∼ 0.38, the value corresponding to the mean stellar mass, to convert our Z into Zg values. In the right panel of Fig. 19 we show the MZRg obtained following the above approach. Despite the fact that our C3 sample is a subsample of the one used by Cullen et al. (2021), the assumed [O/Fe] appears reasonable, as the C3 stacks shown in Fig. 19 follow a consistent trend with the results found by Cullen et al. (2021) using a different stacking procedure.

In order to obtain an independent value for the gas metallicity, we also probed the Si3-O3C3 calibration presented by Byler et al. (2020). For our stacks, the values obtained with the assumed [O/Fe] are consistent within 0.15 dex with the gas-phase metallicities obtained using the Si3-O3C3 calibration. We note that the latter is found to have a median offset of 0.35 dex when compared to other well-known metallicity calibrations based on optical indices. Upon acknowledging these differences and the relatively good agreement between these two methods, we chose to use the re-scaled values with the mean [O/Fe] in the following sections. Clearly, follow up studies probing bright optical lines of CIII] emitters are necessary to provide more reliable gas-phase metallicity determinations.

On the right panel of Fig. 19, our best fit with the assumed mean [O/Fe] is fully consistent with the MZRg at z ∼ 2.3 from Steidel et al. (2014). Our results are also consistent with most data points in the MZRg estimates obtained at z ∼ 3 from Onodera et al. (2016) and Sanders et al. (2021), which are also included for comparison. Our results are thus also consistent with the reported redshift evolution of the MZRg, illustrated here using a comparison with the local relation found by Andrews & Martini (2013). The larger differences are found with respect to the MZRg found by Troncoso et al. (2014), which could be explained by systematic differences in the excitation conditions and metallicities between the samples, as suggested in Sanders et al. (2021). Indeed, the Troncoso et al. (2014) sample is likely to be biased towards lower metallicities. But the agreement of our relation with previous MZRg depends on the α-enhanced assumed in our transformation between metallicity phases.

To better constrain the gas-phase metallicity of CIII] emitters at these redshifts, we need NIR follow-up observations to obtain the rest-optical spectra of these objects. In a future paper, we will present our analysis based on individual galaxies.

Finally, we used the gas-phase metallicity of the C3 sample to probe the Fundamental Metallicity Relation (FMR, Mannucci et al. 2010), which describes an invariant dependence with SFR of the MZRg metallicity of galaxies out to z ∼ 2.5. A recent work by Sanders et al. (2021) suggests that this lack of evolution extends out to z ∼ 3.3. In our work, exploring this relation is highly dependent on the adopted Zg. For example, if we assume that Z = Zg our results would show an offset to low metallicity of ∼0.5 dex. However, assuming the average α-enhancement derived by Cullen et al. (2021) and used in Fig. 19, we find a trend that appears in agreement with the slope of the FMR, as shown in Fig. 20. While two stacks in the B scheme (i.e., only galaxies with z ≥ 2.9), appear offset towards lower metallicity, D stacks (i.e., only galaxies with z < 2.9) and A stacks (stars, representing all galaxies at 2.4 ≲ z ≲ 3.9) find a relatively good agreement with the slope of local FMR. Considering the typically large uncertainties involved both in the data measurements and metallicity derivation, especially those inherent to the different spectral features and methodologies applied in this and previous works, this result is surprisingly robust. Showing itself to be in agreement with recent results (Sanders et al. 2021), it favors the scenario where the FMR does not evolves significantly up to z ∼ 3.

thumbnail Fig. 20.

FMR with the stacks A, B, D by stellar mass, color-coded by EW(CIII]). Blue dashed line is the relation shown in Mannucci et al. (2010). The right scale is the [O/H] ∼ log(Zg/Z) assuming the α-enhanced in Eq. (11).

4.2. Stellar metallicity–C/O relation at z ∼ 3

The C/O ratio may provide us general trends in the evolutionary state of a galaxy and its ISM. In evolved, metal-enriched galaxies, an increase of C/O with increasing metallicity has been observed (Garnett et al. 1995; Berg et al. 2016, 2019a) and also reproduced by models (e.g., Henry et al. 2000; Mollá et al. 2015; Mattsson 2010). This trend can be explained because carbon is primarily produced by the triple-α process in both massive and low- to intermediate-mass stars, however, in massive stars, carbon arises almost exclusively from the production due to metallicity-dependent stellar winds, mass loss, and ISM enrichment – which are greater at higher metallicities (Henry et al. 2000). An evolutionary effect due to the delayed release of carbon (which is mostly produced by low- and intermediate-mass) relative to oxygen (which is produced almost exclusively by massive stars) in younger and less metal-rich systems is an alternative explanation for this trend (Garnett et al. 1995).

In this section, we discuss the position of the C3 sample in the Z–C/O plane, which provide new insights on their chemical enrichment. In the left panel in Fig. 21, we present the Z–C/O relation for our C3 sample in VANDELS along with the predictions drawn from chemical evolution models (model MWG-11 in Romano et al. 2020 as well as model K20 in Kobayashi et al. 2020) and stellar metallicities [Fe/H] estimations of local Milky Way stars in the halo, thick, and thin disk (Amarsi et al. 2019). We show results from stack A, which include the entire C3 sample. Symbols are color-coded by EW(CIII]) and their marker depend on the parameter used for stacking, that is, the stacks based on the EW(CIII]), FUV, and Ks bands, depicted as stars, left-triangle, and pentagon, respectively. Our results indicate that C/O increases with stellar metallicity for Z ≳ 10% solar, in agreement with the trends of chemical evolution model and of metal-rich stars in the Galactic thick disk. Our results show higher values that those predicted by Kobayashi et al. (2020) for a given stellar metallicity, but this also occurs with most of the thin-disk stars. According to these authors, the latter can be partially explained by an under-prediction of carbon yields by AGB stars on the models. More generally, we note that our results are found to be in better agreement with the model prediction by Romano et al. (2020).

thumbnail Fig. 21.

Relation between C/O and metallicity. Left panel: C/O–Z relation. Black symbols are values of stars from the Galactic thin, thick disk, halo and unclassified from the literature (Amarsi et al. 2019). The dashed black line is the K20 model in Kobayashi et al. (2020) and the dotted line is the model MWG-11 in Romano et al. (2020). Right panel: C/O–Zg relation. The red markers are high-redshift galaxies from literature and black markers are local galaxies and HII regions (see text for references). The multi-zone chemical evolution models from Mattsson (2010) are also shown by black lines. The green dotted-dashed line is the chemical model from Mollá et al. (2015). Both panels: the red dashed line is the C/O solar value. The stack A by EW(CIII]) and broad band luminosities are shown by markers according to legend and color-coded by EW(CIII]).

For the gas-phase metallicity, we present results on the right panel of Fig. 21. The stacks are the same on the left panel, but re-scaled based on Eq. (11). We include four different predictions from chemical evolution models by Mattsson (2010) and Mollá et al. (2015). We also include comparison samples from the literature: LAEs (Bayliss et al. 2014; Christensen et al. 2012; Erb et al. 2010; James et al. 2014; Villar-Martín et al. 2004) and LGBs (de Barros et al. 2016; Vanzella et al. 2016; Steidel et al. 2016) at z ∼ 2 − 4, and local analogs such as blue compact dwarfs (BCD; Garnett et al. 1995, 1997; Kobulnicky et al. 1997; Kobulnicky & Skillman 1998; Izotov et al. 1999; Thuan et al. 1999; Berg et al. 2016; Senchyna et al. 2021) and HII regions (Garnett et al. 1995, 1999; Kurt et al. 1995; Mattsson 2010; Senchyna et al. 2021) at z ∼ 0. We note that metallicities for these objects are derived from nebular lines. Our points in Fig. 21, instead, assume Eq. (11) with a fixed [O/Fe]. However, this assumption might not be accurate and significant (and probably different) levels of α-enhancement in individual galaxies may generate systematic offsets between Z and Zg (Cullen et al. 2021). Such effects will tend to increase the dispersion of data in the x-axis. With these caveats in mind, the trends shown in Fig. 21 are still useful with regard to discussing the different levels of chemical enrichment traced by the C/O ratio, particularly for our C3 sample.

Despite the large scatter shown by observations, Fig. 21 shows a trend of increasing C/O with stellar metallicity for our VANDELS C3 sample. This trend is strongly driven by the low mass (low luminosity) stacks, which show the lower C/O and metallicities. Particularly, the stack by EW(CIII]) with the highest EW(CIII]) (purple star in Fig. 21) has the lower stellar metallicity and a lower C/O ratio, suggesting that this population with extreme EW(CIII]) does indeed host a young stellar population still dominated by massive stars.

The comparison of the VANDELS C3 sample with the sample of galaxies at similar redshift (red symbols) shows that our stacks have higher metallicity (even assuming Z = Zg) and higher C/O ratios than that of most CIII] emitters in the comparison sample. This suggests that, on average, we are probing galaxies that are chemically more evolved.

Compared to the local sample, the C/O ratios of our VANDELS C3 sample are comparable to those of local HII regions at similar metallicity, for which the C/O abundances increase toward solar values due to a mix of young and aged stellar populations contributing to C production. Local BCDs, instead, show similar C/O but lower metallicities; that is to say that these are compatible with our results if we assume a lower α-enhancement.

Several mechanisms can affect the position of galaxies in Fig. 21. This includes variations in the star formation histories (i.e., shape, duration of bursts, and star formation efficiency), gas fraction and inflow rates that may affect the gas metallicity and even potential changes in the initial mass function – all of which may contribute to the scatter (e.g., Mattsson 2010; Mollá et al. 2015; Berg et al. 2016; Vincenzo & Kobayashi 2018; Kobayashi et al. 2020; Palla et al. 2020). High values of C/O at low metallicity could be due to the effects of a massive inflow of pristine gas, which may lower the gas metallicity of an otherwise more chemically evolved system, without altering C/O (Nakajima et al. 2018a). This is another caveat to consider when there is not a direct measurement of the gas-phase metallicity.

In models, C/O is sensitive to the prescriptions for yields, the initial mass function, star formation efficiencies, and inflow rates that are not fully understood that make difficult to interpret observations and may also produce variable levels of C/O at low metallicity, as suggested in Mattsson (2010) and Berg et al. (2016). Although models can reproduce some observations, they do not explain the large scatter observed and the precise values completely. Moreover, the large observational and methodological uncertainties involved make the interpretation of the C/O abundances challenging.

In this work, we find that our results are consistent with the trend of increasing C/O seen in the models for the range of gas-phase metallicities probed by the stacks (around 20 to 60% solar), assuming a constant α-enhancement. In particular, we compare our data with the multi-zone chemical evolution models from Mattsson (2010). Case A1 (black dashed line in Fig. 21) represents a case where the low- and intermediate-mass stars are producing most of the carbon, while case B1 (black dotted line in Fig. 21) represents a case where the carbon is to a large extent produced in high mass stars. Other models such as those in Mollá et al. (2015) (see also Berg et al. 2016), which predict higher values of log(C/O) for a given O/H, find a good agreement especially for our stacks at higher EW(CIII]).

Based on a detailed comparison with the models, Berg et al. (2019a) discuss the sensitivity of the C/O ratio to both the detailed star formation history and supernova (SNe) feedback. The models demonstrate that lower C/O values are found at lower star formation efficiencies and burst over a longer duration. On the other hand, larger C/O ratios are more related to the presence of SNe feedback that ejects oxygen and reduces the effective yields. The rapid C/O enrichment found for the C3 sample could be thus related to a more bursty star formation history. Given that our C3 sample is mostly made up of main-sequence galaxies, we speculate that a recent star formation history with multiple bursts could have enhanced the C/O to their observed levels.

More detailed observations are needed for both local and high redshift galaxies to fully understand the C/O enrichment history through cosmic time. The HST COS Legacy Archive Spectroscopic SurveY (CLASSY; Berg et al. 2019b), which has obtained high S/N UV nebular spectra for a representative sample of local galaxies, will certainly contribute to clarifying this scenario in the local universe. In addition, future observations with JWST will enable more detailed studies of the C/O evolution of star-forming galaxies at intermediate and high redshift.

5. Conclusions

In this work, we studied a broad representative sample of 217 galaxies with CIII] detection at 2 < z < 4, covering a range of ∼2 dex in stellar mass. These CIII] emitters have a broad range of UV luminosities, thus allowing for a stacking analysis to characterize their mean physical properties. We considered stacking by different bins of stellar mass, rest-frame FUV, and Ks luminosities, as well as rest-frame EW(CIII]) and EW(Lyα). We derived the stellar metallicity and the C/O abundance from stack spectra. In this paper, we discuss several relations found between these parameters. We summarize our conclusions as follows:

  • (i)

    Reliable (S/N > 3) CIII] emitters represent ∼30% of the VANDELS parent sample at z ∼ 2 − 4. They show EW(CIII]) between 0.3 and 20 Å with a mean value of ∼4 Å. However, stacked spectra of galaxies with marginal or non detection of CIII] (S/N < 3) in individual spectra show weak (EW ≲ 2 Å) CIII] emission, suggesting this line is common in normal star-forming galaxies at z ∼ 3. On the other hand, extreme emitters (EW(CIII]) ≳ 8 Å) are exceedingly rare (∼3%) in VANDELS, which is expected as the C3 sample is drawn from a parent sample of main-sequence galaxies.

  • (ii)

    Stacking reveals that faint UV nebular lines of OIII], SiIII, CIV, and HeII, as well as fine-structure emission lines of SiII are ubiquitous in CIII] emitters. We find that the strength of the nebular lines depends on the stellar mass and luminosity. Overall, less massive (and fainter in any band) stacks show more intense Lyα, CIV, HeII, OIII] and CIII] than more massive (and brighter in any band), which tend to have fainter UV emission lines.

  • (iii)

    According to UV diagnostic diagrams, nebular lines in VANDELS CIII] emitters are powered by stellar photoionization, suggesting no ionization source other than massive stars. These results add new constraints to models on the flux ratios and EWs of nebular emission lines in main-sequence SF regions for the faint CIII] regime at z ∼ 2 − 4.

  • (iv)

    The EW(Lyα) and EW(CIII]) appear correlated in our C3 sample but a large scatter is found for individual galaxies. Stacks with larger EW(CIII]) show larger EW(Lyα), but not all CIII] emitters are Lyα emitters.

  • (v)

    Galaxies with higher EW(CIII]) show lower stellar metallicities. This result suggests that extremely low metallicities (< 10% solar) should be expected for the most extreme galaxies in terms of their EW(CIII]). Our results (Eq. (3)) show that galaxies with stellar metallicity Z < 0.1 Z are typically strong CIII] emitters (EW(CIII]) > 7 Å).

  • (vi)

    The stellar metallicities of CIII] emitters are not significantly different from that of the parent sample, increasing from ∼10% to ∼40% solar for stellar masses log(M/M) ∼ 9−10.5. The stellar mass-metallicity relation of the CIII] emitters is consistent with previous works showing strong evolution from z = 0 to z ∼ 3, both using stellar (Fe/H) and nebular (O/H) metallicities inferred from the UV spectra and assuming an average [Fe/O] = 0.38, recently found for a subsample of VANDELS galaxies at z ∼ 3 (Cullen et al. 2021).

  • (vii)

    We find the C/O abundances of CIII] emitters ranging between 35%−150% solar, with a noticeable increase with FUV luminosity, and a smooth decrease with CIII] EWs. Fainter FUV galaxies have lower C/O, higher EW(CIII]), and lower Z, which suggest a UV spectrum dominated by massive stars and a bright nebular component which is still chemically unevolved.

  • (viii)

    We discuss, for the first time, the C/O–Fe/H and the C/O–O/H relations for star-forming galaxies at z ∼ 3. They show stellar and nebular abundances consistent with the trends observed in Milky Way halo and thick-disk stars and local HII galaxies, respectively. We find a good agreement with modern chemical evolution models, which suggest that CIII] emitters at z ∼ 3 are experiencing an active phase of chemical enrichment.

Our results provide new insights into the nature of UV line emitters at z ∼ 3, paving the way for future studies at higher-z using the James Webb Space Telescope.


4

Using the extinction code at http://github.com/kbarbary/extinction

5

We used the version 3.2 publicly available at https://www.iaa.csic.es/~epm/HII-CHI-mistry-UV.html

Acknowledgments

We thank the anonymous referee for the detailed review and useful suggestions. This work is based on data products from observations made with ESO Telescopes at La Silla Paranal Observatory under ESO program ID 194.A-2003 (PIs: Laura Pentericci and Ross McLure). MLl acknowledges support from the National Agency for Research and Development (ANID)/Scholarship Program/Doctorado Nacional/2019-21191036. RA acknowledges support from ANID FONDECYT Regular Grant 1202007. This work has made extensive use of Python packages astropy (Astropy Collaboration 2018), numpy (Harris et al. 2020), and Matplotlib (Hunter 2007).

References

  1. Amarsi, A. M., Nissen, P. E., & Skúladóttir, Á. 2019, A&A, 630, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Amorín, R. O., Pérez-Montero, E., & Vílchez, J. M. 2010, ApJ, 715, L128 [CrossRef] [Google Scholar]
  3. Amorín, R., Pérez-Montero, E., Vílchez, J. M., & Papaderos, P. 2012, ApJ, 749, 185 [CrossRef] [Google Scholar]
  4. Amorín, R., Pérez-Montero, E., Contini, T., et al. 2015, A&A, 578, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Amorín, R., Fontana, A., Pérez-Montero, E., et al. 2017, Nat. Astron., 1, 0052 [Google Scholar]
  6. Andrews, B. H., & Martini, P. 2013, ApJ, 765, 140 [NASA ADS] [CrossRef] [Google Scholar]
  7. Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009, ARA&A, 47, 481 [NASA ADS] [CrossRef] [Google Scholar]
  8. Astropy Collaboration (Price-Whelan, A. M., et al.) 2018, AJ, 156, 123 [Google Scholar]
  9. Bacon, R., Conseil, S., Mary, D., et al. 2017, A&A, 608, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Baldwin, J. A., Phillips, M. M., & Terlevich, R. 1981, PASP, 93, 5 [Google Scholar]
  11. Bayliss, M. B., Rigby, J. R., Sharon, K., et al. 2014, ApJ, 790, 144 [NASA ADS] [CrossRef] [Google Scholar]
  12. Berg, D. A., Skillman, E. D., Henry, R. B. C., Erb, D. K., & Carigi, L. 2016, ApJ, 827, 126 [NASA ADS] [CrossRef] [Google Scholar]
  13. Berg, D. A., Erb, D. K., Auger, M. W., Pettini, M., & Brammer, G. B. 2018, ApJ, 859, 164 [CrossRef] [Google Scholar]
  14. Berg, D. A., Erb, D. K., Henry, R. B. C., Skillman, E. D., & McQuinn, K. B. W. 2019a, ApJ, 874, 93 [NASA ADS] [CrossRef] [Google Scholar]
  15. Berg, D., Chisholm, J., Heckman, T. M., et al. 2019b, The COS Legacy Archive Spectroscopic SurveY (CLASSY): A UV Treasury of Star-Forming Galaxies, 15840 [Google Scholar]
  16. Bouwens, R. J., Smit, R., Labbé, I., et al. 2016, ApJ, 831, 176 [Google Scholar]
  17. Bressan, A., Marigo, P., Girardi, L., et al. 2012, MNRAS, 427, 127 [Google Scholar]
  18. Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000 [NASA ADS] [CrossRef] [Google Scholar]
  19. Byler, N., Dalcanton, J. J., Conroy, C., et al. 2018, ApJ, 863, 14 [NASA ADS] [CrossRef] [Google Scholar]
  20. Byler, N., Kewley, L. J., Rigby, J. R., et al. 2020, ApJ, 893, 1 [NASA ADS] [CrossRef] [Google Scholar]
  21. Calabrò, A., Castellano, M., Pentericci, L., et al. 2021, A&A, 646, A39 [EDP Sciences] [Google Scholar]
  22. Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682 [NASA ADS] [CrossRef] [Google Scholar]
  23. Cardamone, C., Schawinski, K., Sarzi, M., et al. 2009, MNRAS, 399, 1191 [NASA ADS] [CrossRef] [Google Scholar]
  24. Carigi, L. 1994, ApJ, 424, 181 [NASA ADS] [CrossRef] [Google Scholar]
  25. Carigi, L., & Peimbert, M. 2011, Rev. Mex. Astron. Astrofis., 47, 139 [NASA ADS] [Google Scholar]
  26. Carnall, A. C. 2017, ArXiv e-prints [arXiv:1705.05165] [Google Scholar]
  27. Carnall, A. C., McLure, R. J., Dunlop, J. S., & Davé, R. 2018, MNRAS, 480, 4379 [Google Scholar]
  28. Cassata, P., Tasca, L. A. M., Le Fèvre, O., et al. 2015, A&A, 573, A24 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  29. Castellano, M., Sommariva, V., Fontana, A., et al. 2014, A&A, 566, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  30. Chabrier, G. 2003, PASP, 115, 763 [Google Scholar]
  31. Chiappini, C., Romano, D., & Matteucci, F. 2003, MNRAS, 339, 63 [NASA ADS] [CrossRef] [Google Scholar]
  32. Christensen, L., Laursen, P., Richard, J., et al. 2012, MNRAS, 427, 1973 [NASA ADS] [CrossRef] [Google Scholar]
  33. Cullen, F., McLure, R. J., Dunlop, J. S., et al. 2019, MNRAS, 487, 2038 [Google Scholar]
  34. Cullen, F., McLure, R. J., Dunlop, J. S., et al. 2020, MNRAS, 495, 1501 [Google Scholar]
  35. Cullen, F., Shapley, A. E., McLure, R. J., et al. 2021, MNRAS, 505, 903 [CrossRef] [Google Scholar]
  36. Dayal, P., Volonteri, M., Choudhury, T. R., et al. 2020, MNRAS, 495, 3065 [Google Scholar]
  37. de Barros, S., Vanzella, E., Amorín, R., et al. 2016, A&A, 585, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Du, X., Shapley, A. E., Tang, M., et al. 2020, ApJ, 890, 65 [NASA ADS] [CrossRef] [Google Scholar]
  39. Endsley, R., Stark, D. P., Chevallard, J., & Charlot, S. 2021, MNRAS, 500, 5229 [Google Scholar]
  40. Erb, D. K., Shapley, A. E., Pettini, M., et al. 2006, ApJ, 644, 813 [NASA ADS] [CrossRef] [Google Scholar]
  41. Erb, D. K., Pettini, M., Shapley, A. E., et al. 2010, ApJ, 719, 1168 [Google Scholar]
  42. Esteban, C., García-Rojas, J., Carigi, L., et al. 2014, MNRAS, 443, 624 [NASA ADS] [CrossRef] [Google Scholar]
  43. Falcón-Barroso, J., Sánchez-Blázquez, P., Vazdekis, A., et al. 2011, A&A, 532, A95 [Google Scholar]
  44. Fan, X., Strauss, M. A., Becker, R. H., et al. 2006, AJ, 132, 117 [NASA ADS] [CrossRef] [Google Scholar]
  45. Feltre, A., Charlot, S., & Gutkin, J. 2016, MNRAS, 456, 3354 [Google Scholar]
  46. Feltre, A., Maseda, M. V., Bacon, R., et al. 2020, A&A, 641, A118 [EDP Sciences] [Google Scholar]
  47. Finkelstein, S. L., D’Aloisio, A., Paardekooper, J.-P., et al. 2019, ApJ, 879, 36 [Google Scholar]
  48. Fuller, S., Lemaux, B. C., Bradač, M., et al. 2020, ApJ, 896, 156 [NASA ADS] [CrossRef] [Google Scholar]
  49. Gallazzi, A., Charlot, S., Brinchmann, J., White, S. D. M., & Tremonti, C. A. 2005, MNRAS, 362, 41 [Google Scholar]
  50. Garilli, B., McLure, R., Pentericci, L., et al. 2021, A&A, 647, A150 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Garnett, D. R., Skillman, E. D., Dufour, R. J., et al. 1995, ApJ, 443, 64 [Google Scholar]
  52. Garnett, D. R., Skillman, E. D., Dufour, R. J., & Shields, G. A. 1997, ApJ, 481, 174 [CrossRef] [Google Scholar]
  53. Garnett, D. R., Shields, G. A., Peimbert, M., et al. 1999, ApJ, 513, 168 [NASA ADS] [CrossRef] [Google Scholar]
  54. Grazian, A., Giallongo, E., Boutsia, K., et al. 2018, A&A, 613, A44 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  55. Guseva, N. G., Izotov, Y. I., Schaerer, D., et al. 2020, MNRAS, 497, 4293 [CrossRef] [Google Scholar]
  56. Gutkin, J., Charlot, S., & Bruzual, G. 2016, MNRAS, 462, 1757 [Google Scholar]
  57. Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [Google Scholar]
  58. Henry, R. B. C., Edmunds, M. G., & Köppen, J. 2000, ApJ, 541, 660 [NASA ADS] [CrossRef] [Google Scholar]
  59. Hirschmann, M., Charlot, S., Feltre, A., et al. 2019, MNRAS, 487, 333 [NASA ADS] [CrossRef] [Google Scholar]
  60. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [Google Scholar]
  61. Hutchison, T. A., Papovich, C., Finkelstein, S. L., et al. 2019, ApJ, 879, 70 [NASA ADS] [CrossRef] [Google Scholar]
  62. Izotov, Y. I., Chaffee, F. H., Foltz, C. B., et al. 1999, ApJ, 527, 757 [NASA ADS] [CrossRef] [Google Scholar]
  63. Izotov, Y. I., Schaerer, D., Thuan, T. X., et al. 2016, MNRAS, 461, 3683 [Google Scholar]
  64. James, B. L., Pettini, M., Christensen, L., et al. 2014, MNRAS, 440, 1794 [NASA ADS] [CrossRef] [Google Scholar]
  65. Jaskot, A. E., & Ravindranath, S. 2016, ApJ, 833, 136 [NASA ADS] [CrossRef] [Google Scholar]
  66. Kewley, L. J., Nicholls, D. C., & Sutherland, R. S. 2019, ARA&A, 57, 511 [Google Scholar]
  67. Kobayashi, C., Karakas, A. I., & Lugaro, M. 2020, ApJ, 900, 179 [Google Scholar]
  68. Kobulnicky, H. A., & Skillman, E. D. 1998, ApJ, 497, 601 [NASA ADS] [CrossRef] [Google Scholar]
  69. Kobulnicky, H. A., Skillman, E. D., Roy, J.-R., Walsh, J. R., & Rosa, M. R. 1997, ApJ, 477, 679 [NASA ADS] [CrossRef] [Google Scholar]
  70. Kocevski, D. D., Hasinger, G., Brightman, M., et al. 2018, ApJS, 236, 48 [CrossRef] [Google Scholar]
  71. Kurt, C. M., Dufour, R. J., Garnett, D. R., et al. 1995, in Rev. Mex. Astron. Astrofis. Conf. Ser., eds. M. Pena, & S. Kurtz, 3, 223 [NASA ADS] [Google Scholar]
  72. Laporte, N., Nakajima, K., Ellis, R. S., et al. 2017, ApJ, 851, 40 [NASA ADS] [CrossRef] [Google Scholar]
  73. Le Fèvre, O., Tasca, L. A. M., Cassata, P., et al. 2015, A&A, 576, A79 [Google Scholar]
  74. Le Fèvre, O., Lemaux, B. C., Nakajima, K., et al. 2019, A&A, 625, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Leitherer, C., Tremonti, C. A., Heckman, T. M., & Calzetti, D. 2011, AJ, 141, 37 [NASA ADS] [CrossRef] [Google Scholar]
  76. Leitherer, C., Ekström, S., Meynet, G., et al. 2014, ApJS, 212, 14 [Google Scholar]
  77. Lian, J. H., Li, J. R., Yan, W., & Kong, X. 2015, MNRAS, 446, 1449 [NASA ADS] [CrossRef] [Google Scholar]
  78. Luo, B., Brandt, W. N., Xue, Y. Q., et al. 2017, ApJS, 228, 2 [Google Scholar]
  79. Mainali, R., Kollmeier, J. A., Stark, D. P., et al. 2017, ApJ, 836, L14 [Google Scholar]
  80. Maiolino, R., & Mannucci, F. 2019, A&ARv, 27, 3 [Google Scholar]
  81. Maiolino, R., Nagao, T., Grazian, A., et al. 2008, A&A, 488, 463 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  82. Mannucci, F., Cresci, G., Maiolino, R., Marconi, A., & Gnerucci, A. 2010, MNRAS, 408, 2115 [NASA ADS] [CrossRef] [Google Scholar]
  83. Marchi, F., Pentericci, L., Guaita, L., et al. 2017, A&A, 601, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  84. Marchi, F., Pentericci, L., Guaita, L., et al. 2019, A&A, 631, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  85. Marigo, P., Bressan, A., Nanni, A., Girardi, L., & Pumo, M. L. 2013, MNRAS, 434, 488 [Google Scholar]
  86. Maseda, M. V., van der Wel, A., Rix, H.-W., et al. 2014, ApJ, 791, 17 [CrossRef] [Google Scholar]
  87. Maseda, M. V., Brinchmann, J., Franx, M., et al. 2017, A&A, 608, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  88. Mason, C. A., Fontana, A., Treu, T., et al. 2019, MNRAS, 485, 3947 [NASA ADS] [CrossRef] [Google Scholar]
  89. Matthee, J., Sobral, D., Hayes, M., et al. 2021, MNRAS, 505, 1382 [NASA ADS] [CrossRef] [Google Scholar]
  90. Mattsson, L. 2010, A&A, 515, A68 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  91. McLure, R. J., Pentericci, L., Cimatti, A., et al. 2018, MNRAS, 479, 25 [NASA ADS] [Google Scholar]
  92. Mollá, M., García-Vargas, M. L., & Bressan, A. 2009, MNRAS, 398, 451 [CrossRef] [Google Scholar]
  93. Mollá, M., Cavichia, O., Gavilán, M., & Gibson, B. K. 2015, MNRAS, 451, 3693 [CrossRef] [Google Scholar]
  94. Naidu, R. P., Tacchella, S., Mason, C. A., et al. 2020, ApJ, 892, 109 [NASA ADS] [CrossRef] [Google Scholar]
  95. Nakajima, K., Schaerer, D., Le Fèvre, O., et al. 2018a, A&A, 612, A94 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  96. Nakajima, K., Fletcher, T., Ellis, R. S., Robertson, B. E., & Iwata, I. 2018b, MNRAS, 477, 2098 [NASA ADS] [CrossRef] [Google Scholar]
  97. Newville, M., Stensitzki, T., Allen, D. B., & Ingargiola, A. 2014, https://doi.org/10.5281/zenodo.11813 [Google Scholar]
  98. Onodera, M., Carollo, C. M., Lilly, S., et al. 2016, ApJ, 822, 42 [NASA ADS] [CrossRef] [Google Scholar]
  99. Ouchi, M., Ono, Y., & Shibuya, T. 2020, ARA&A, 58, 617 [Google Scholar]
  100. Palla, M., Calura, F., Matteucci, F., et al. 2020, MNRAS, 494, 2355 [NASA ADS] [CrossRef] [Google Scholar]
  101. Paoletti, D., Hazra, D. K., Finelli, F., & Smoot, G. F. 2020, JCAP, 2020, 005 [Google Scholar]
  102. Peña-Guerrero, M. A., Leitherer, C., de Mink, S., Wofford, A., & Kewley, L. 2017, ApJ, 847, 107 [CrossRef] [Google Scholar]
  103. Pentericci, L., Vanzella, E., Fontana, A., et al. 2014, ApJ, 793, 113 [NASA ADS] [CrossRef] [Google Scholar]
  104. Pentericci, L., McLure, R. J., Garilli, B., et al. 2018, A&A, 616, A174 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  105. Pérez-Montero, E., & Amorín, R. 2017, MNRAS, 467, 1287 [Google Scholar]
  106. Pérez-Montero, E., Contini, T., Lamareille, F., et al. 2013, A&A, 549, A25 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  107. Pettini, M., Steidel, C. C., Adelberger, K. L., Dickinson, M., & Giavalisco, M. 2000, ApJ, 528, 96 [Google Scholar]
  108. Ravindranath, S., Monroe, T., Jaskot, A., Ferguson, H. C., & Tumlinson, J. 2020, ApJ, 896, 170 [NASA ADS] [CrossRef] [Google Scholar]
  109. Reddy, N. A., Kriek, M., Shapley, A. E., et al. 2015, ApJ, 806, 259 [NASA ADS] [CrossRef] [Google Scholar]
  110. Rigby, J. R., Bayliss, M. B., Gladders, M. D., et al. 2015, ApJ, 814, L6 [NASA ADS] [CrossRef] [Google Scholar]
  111. Robertson, B. E., Ellis, R. S., Furlanetto, S. R., & Dunlop, J. S. 2015, ApJ, 802, L19 [Google Scholar]
  112. Romano, D., Franchini, M., Grisoni, V., et al. 2020, A&A, 639, A37 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  113. Runco, J. N., Shapley, A. E., Sanders, R. L., et al. 2021, MNRAS, 502, 2600 [NASA ADS] [CrossRef] [Google Scholar]
  114. Salim, S., Boquien, M., & Lee, J. C. 2018, ApJ, 859, 11 [Google Scholar]
  115. Sanders, R. L., Shapley, A. E., Jones, T., et al. 2021, ApJ, 914, 19 [CrossRef] [Google Scholar]
  116. Saxena, A., Pentericci, L., Mirabelli, M., et al. 2020, A&A, 636, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  117. Schaerer, D., Izotov, Y. I., Nakajima, K., et al. 2018, A&A, 616, L14 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  118. Schmidt, K. B., Treu, T., Bradač, M., et al. 2016, ApJ, 818, 38 [NASA ADS] [CrossRef] [Google Scholar]
  119. Senchyna, P., Stark, D. P., Vidal-García, A., et al. 2017, MNRAS, 472, 2608 [NASA ADS] [CrossRef] [Google Scholar]
  120. Senchyna, P., Stark, D. P., Charlot, S., et al. 2021, MNRAS, 503, 6112 [NASA ADS] [CrossRef] [Google Scholar]
  121. Shapley, A. E. 2011, ARA&A, 49, 525 [NASA ADS] [CrossRef] [Google Scholar]
  122. Shapley, A. E., Steidel, C. C., Pettini, M., & Adelberger, K. L. 2003, ApJ, 588, 65 [Google Scholar]
  123. Sobral, D., Matthee, J., Darvish, B., et al. 2015, ApJ, 808, 139 [Google Scholar]
  124. Sommariva, V., Mannucci, F., Cresci, G., et al. 2012, A&A, 539, A136 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  125. Speagle, J. S., Steinhardt, C. L., Capak, P. L., & Silverman, J. D. 2014, ApJS, 214, 15 [Google Scholar]
  126. Stark, D. P. 2016, ARA&A, 54, 761 [Google Scholar]
  127. Stark, D. P., Richard, J., Siana, B., et al. 2014, MNRAS, 445, 3200 [NASA ADS] [CrossRef] [Google Scholar]
  128. Stark, D. P., Richard, J., Charlot, S., et al. 2015, MNRAS, 450, 1846 [NASA ADS] [CrossRef] [Google Scholar]
  129. Stark, D. P., Ellis, R. S., Charlot, S., et al. 2017, MNRAS, 464, 469 [NASA ADS] [CrossRef] [Google Scholar]
  130. Steidel, C. C., Giavalisco, M., Pettini, M., Dickinson, M., & Adelberger, K. L. 1996, ApJ, 462, L17 [Google Scholar]
  131. Steidel, C. C., Pettini, M., & Adelberger, K. L. 2001, ApJ, 546, 665 [Google Scholar]
  132. Steidel, C. C., Rudie, G. C., Strom, A. L., et al. 2014, ApJ, 795, 165 [Google Scholar]
  133. Steidel, C. C., Strom, A. L., Pettini, M., et al. 2016, ApJ, 826, 159 [NASA ADS] [CrossRef] [Google Scholar]
  134. Tang, M., Stark, D. P., Chevallard, J., et al. 2021, MNRAS, 501, 3238 [NASA ADS] [CrossRef] [Google Scholar]
  135. Thuan, T. X., Izotov, Y. I., & Foltz, C. B. 1999, ApJ, 525, 105 [NASA ADS] [CrossRef] [Google Scholar]
  136. Topping, M. W., Shapley, A. E., Reddy, N. A., et al. 2020, MNRAS, 495, 4430 [Google Scholar]
  137. Tremonti, C. A., Heckman, T. M., Kauffmann, G., et al. 2004, ApJ, 613, 898 [Google Scholar]
  138. Troncoso, P., Maiolino, R., Sommariva, V., et al. 2014, A&A, 563, A58 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  139. Vanzella, E., De Barros, S., Cupani, G., et al. 2016, ApJ, 821, L27 [Google Scholar]
  140. Vanzella, E., Caminha, G. B., Rosati, P., et al. 2021, A&A, 646, A57 [EDP Sciences] [Google Scholar]
  141. Villar-Martín, M., Cerviño, M., & González Delgado, R. M. 2004, MNRAS, 355, 1132 [CrossRef] [Google Scholar]
  142. Vincenzo, F., & Kobayashi, C. 2018, A&A, 610, L16 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  143. Wang, B., Heckman, T. M., Amorín, R., et al. 2021, ApJ, 916, 3 [NASA ADS] [CrossRef] [Google Scholar]
  144. Wise, J. H., Demchenko, V. G., Halicek, M. T., et al. 2014, MNRAS, 442, 2560 [NASA ADS] [CrossRef] [Google Scholar]
  145. Yang, J., Wang, F., Fan, X., et al. 2020, ApJ, 904, 26 [Google Scholar]
  146. Yue, B., Castellano, M., Ferrara, A., et al. 2018, ApJ, 868, 115 [NASA ADS] [CrossRef] [Google Scholar]
  147. Zahid, H. J., Kudritzki, R.-P., Conroy, C., Andrews, B., & Ho, I. T. 2017, ApJ, 847, 18 [Google Scholar]

Appendix A: Analysis of non-detected CIII] emitters in the parent sample

The results presented in this paper are based on the selection of a sample of galaxies with the CIII] line detected with S/N> 3. This S/N threshold is motivated by the need of having reliable systematic redshifts, which are obtained fitting this nebular line –the most intense UV metal line in the parent sample. In this section, we study whether this selection criterion has any effect on our results. We analyze the 521 star-forming galaxies from the VANDELS parent sample with non-detections in the CIII] line, namely, those objects excluded from our C3 sample. We performed the same stacking analysis done for the C3 sample, following schemes A and B in bins of FUV luminosity and using the same bin distribution and methodology explained in Section 2. Because of the lack of accurate systemic redshifts for this sample, stacks were performed using the spectroscopic redshift reported by the VANDELS DR3 catalog, which is mostly based on template fitting of ISM absorption features by the tool PANDORA/EZ (McLure et al. 2018). We find that the resulting stacks show similar features compared to the stacks in the C3 sample, showing detections of multiple ISM absorption lines and nebular emission lines. Overall, the effect of using spectroscopic instead of systemic redshifts in the stacking is to broaden emission line profiles in the stacks due to small velocity shifts. As a result, the S/N ratio of the detected lines tend to be lower and their line fluxes more uncertain than that of the C3 sample, despite the fact that the number of galaxies averaged in each stack is generally higher.

In all the stacks but the one at lower luminosity, we detect CIII] with a S/N> 3. However, their fluxes are fainter than their counterparts in the C3 sample. In Fig. A.1 we present the resulting EW(CIII]) as a function of the mean FUV luminosity for the stacks A and B with detected CIII]. Their trend is in good agreement with the same relations found using the C3 sample, showing a very mild increase of EW(CIII]) from ∼1Å to ∼2Å towards lower FUV luminosity. This suggests that including these galaxies without individual detection of CIII] in the stacks for the C3 sample will not change the resulting EWs and their relations with luminosity. A similar conclusion is found when comparing line ratios of UV metal lines. While not all the stacks in this test have good detections (S/N> 3) in all the relevant emission lines, the position of their line ratios in the diagnostic diagrams presented in Fig. 11 and 12 are consistent with those of the faintest CIII] emitters, that is, they are fully consistent with pure stellar photoionization.

thumbnail Fig. A.1.

MFUV–EW(CIII]) relation for the stacks by FUV luminosity with the non-detected CIII] sample (in cyan symbols). Faint symbols are the same as in Fig. 10.

Finally, we also test the possible effects in the C/O abundances. Following the same methods described in Section 3.6 for the C3 sample, we estimate C/O in the non-C3 stacks where we have at least S/N≥ 2 detections of OIII]. We find consistent C/O values with those of the C3 sample, with log(C/O) ranging from −0.6 to −0.34 (i.e., 45% to 80% solar, respectively). These values show good agreement with the relation between C/O and EW(CIII]) found in Eq. (5) for the C3 sample, as illustrated in Fig. A.2. Again, this result suggests that the non-detected CIII] emitters in the parent sample that were excluded from the stacking analysis of the C3 sample include intrinsically faint CIII] emitters, which are well represented by the results found for the C3 sample throughout this paper. We can therefore conclude that CIII] emission should be a common feature in UV faint galaxies at z∼ 3. This result is relevant for reionization studies because the number of UV faint galaxies increase significantly at higher redshift (e.g., Bouwens et al. 2016; Yue et al. 2018).

thumbnail Fig. A.2.

EW(CIII])–C/O relation for the stacks by FUV luminosity with the non-detected CIII] sample (in cyan symbols). Faint symbols, shaded regions, and dashed lines are the same as in Fig. 16

Appendix B: Stacked spectra

In this section, we present the resulting stacks B (Fig. B.1), C (Fig. B.2), and D (Fig. B.3) grouping galaxies by FUV and Ks luminosities, and by stellar mass, as well as stacks by EW(CIII]) (Fig. B.4), which are used in the different analyses of this work.

thumbnail Fig. B.1.

Resulting Stack B for each physical parameter. Left to right: MKs, MFUV, and stellar mass. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line in the 1-σ error spectrum. The vertical lines mark known UV lines. Information about the number of galaxies, the mean redshift, and the mean parameter are included in each panel.

thumbnail Fig. B.2.

Resulting Stack C for each physical parameter. Left to right: MKs, MFUV, and stellar mass. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line in the 1-σ error spectrum. The vertical lines mark known UV lines. Information about the number of galaxies, the mean redshift and the mean parameter are included in each panel.

thumbnail Fig. B.3.

Resulting Stack D for each physical parameter. Left to right: MKs, MFUV, and stellar mass. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line in the 1-σ error spectrum. The vertical lines mark known UV lines. Information about the number of galaxies, the mean redshift and the mean parameter are included in each panel.

thumbnail Fig. B.4.

Stacks by EW(CIII]). Left to the right on top row: Stack A, Stack D. Middle row: Stack B. Bottom row: Stack C. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line in the 1-σ error spectrum. The vertical lines mark known UV lines. Information about the number of galaxies, the mean redshift and the mean parameter are included in each panel.

Appendix C: Caveats related with stacking methods

In this section, we discuss possible effects of stacking and the binning choice in the relations shown in Fig. 10. For this aim, in Fig. C.1 we compare the EW(CIII]) measured on the stacked spectra with the median EW(CIII]) of the galaxies in each bin as measured in individual spectra. We do this comparison with stacks A, which are obtained after binning by stellar mass and FUV luminosity. The error bars for the median EWs represent the standard deviation in the bin. We notice that the trends are consistent between both methods, with some differences in the actual value of EW(CIII]), but consistent within the large errors. Therefore, measurements of the stacks are representative of the median properties of the galaxies in each bin and we do not expect strong biases due to bin size or the particular choose of the limits.

thumbnail Fig. C.1.

Relation between EW(CIII]) and stellar mass (left panel) and FUV (right panel) luminosity, similar as shown in Fig. 10 for stack A (blue stars) and C3 sample (small blue circles). The cyan pentagons are the median values for the individual galaxies in the bin and their errorbars correspond to the standard deviation in the bin (boundaries marked by vertical dashed line).

On the other hand, we note that the parameter used for stacking may result in different ranges in stellar metallicity and C/O (e.g., Figs. 14, 15, 16, and 18). This is the result of the distribution of the galaxies in the bins depending on their physical properties. As shown in Fig. 3, the relations between stellar mass and luminosities show scatter and then the same galaxies are not in the corresponding same bin for each parameter. However, the trends are still found if we use only one set of stacking. For instance, in Fig. C.2, we use the cyan dashed line to show the relations found if we only considered the stacks by luminosity. We note that the cyan lines follow the same trend that the black lines, which are the relations reported in the paper including the other set of stacks in the legends. They also are within the 3-σ uncertainty of the reported relations. In the EW(Lyα)–Z relation (top-right panel in Fig. C.2), we note that the dynamical range in EW(Lyα) is shorter compared with the entire set of stacks because only upper limits in stellar metallicity are determined due to the lower S/N of the stacks with higher EWs in this set of stacks. The inclusion of additional sets of stacks binned by different parameters allows exploration of a larger dynamical space, generally at the cost of a larger variance, but without affecting the trends found in the relations.

thumbnail Fig. C.2.

Changes in linear fit in relations presented in Figs. 14, 15, 16, and 18, considering only stacks based on luminosities. Symbols are the same as in Figs. 14, 15, 16, and 18, respectively. The cyan dashed line shows the best fit including only the A stacks binned by luminosity. The stacks excluded in the fit are shown using transparent symbols.

Appendix D: Additional tables

In this section, we present Tables D.1D.5 with the corresponding properties of each stack by stellar mass, Ks, and FUV luminosities, EW(Lyα), and EW(CIII]), respectively. In particular, the observed fluxes measured for the lines used in the analysis, and the C/O and Z estimated are reported.

Table D.1.

Parameters estimated for the stacks by stellar mass

Table D.2.

Parameters estimated for the stacks by Ks-band luminosity

Table D.3.

Parameters estimated for the stacks by FUV luminosity

Table D.4.

Parameters estimated for the stacks by EW(Lyα)

Table D.5.

Parameters estimated for the stacks by EW(CIII])

All Tables

Table 1.

Bin ranges used for the stack analysis with BAGPIPES parameters.

Table 2.

Bin ranges used for the stack analysis with EW.

Table D.1.

Parameters estimated for the stacks by stellar mass

Table D.2.

Parameters estimated for the stacks by Ks-band luminosity

Table D.3.

Parameters estimated for the stacks by FUV luminosity

Table D.4.

Parameters estimated for the stacks by EW(Lyα)

Table D.5.

Parameters estimated for the stacks by EW(CIII])

All Figures

thumbnail Fig. 1.

Spectrum of UDS20394, one of the more intense CIII] emitters in the C3 sample, whose estimated parameters are log(M/M) = 9.32, SFR = 4.01 M yr−1, MFUV = −20, MKs = −20.46, EW(Lyα) = 19.2 Å, and EW(CIII]) = 12.3 Å. The green faint line is the de-redshifted VANDELS spectrum and the black line is the same but resampled by a factor of 2. The blue line in the upper panels is the error spectrum. The red line in the intermediate panels is the scaled sky spectrum.

In the text
thumbnail Fig. 2.

Systemic redshift distribution of the CIII] emitting galaxies in the C3 sample.

In the text
thumbnail Fig. 3.

Relations between the resulting BAGPIPES parameters in the VANDELS main-sequence galaxies. The galaxies of the C3 sample are color-coded by star formation rate. Gray points are VANDELS galaxies of the parent sample that were not included in the C3 sample.

In the text
thumbnail Fig. 4.

M–SFR relation with our sample color-coded by EW (CIII]). The orange solid line is the main sequence at z = 3, according to Speagle et al. (2014). Gray points are as in Fig. 3.

In the text
thumbnail Fig. 5.

Distribution of the resulting BAGPIPES parameters in the C3 sample. Left to right: stellar mass, FUV, and Ks band luminosity. The vertical dashed lines represent the ranges of each bin for stacking according to Table 1.

In the text
thumbnail Fig. 6.

EW(CIII]) (left) and EW(Lyα) (right) distributions of the galaxies in the C3 sample. The vertical dashed lines are the limits of the bins used for stacking according to Table 2.

In the text
thumbnail Fig. 7.

Resulting stacked spectra for each physical parameter with the C3 sample (what we call Stack A, see text for more details). Left to right: MKs, MFUV, and stellar mass. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line is the 1-σ error spectrum. The vertical lines mark known UV lines (in black: emission lines; in red: ISM absorption lines). Information about the number of galaxies, the mean redshift, and the mean parameter are included in each panel.

In the text
thumbnail Fig. 8.

Resulting stacks by EW(Lyα). Details are the same as in Fig. 7. Information about the number of galaxies, the mean redshift, and the mean EW(Lyα) are included in each panel. Left panels: relative Lyα strength in each bin.

In the text
thumbnail Fig. 9.

Emission-line flux and EW measurements in the stack by EW(Lyα) > 20. Panel a: continuum-subtracted spectrum. The orange line shows the error spectrum (1σ). The blue shaded region are the pixel integrated for the CIV line. Panels b and c: red line is the slinefit fitting for the spectrum. The orange shaded region shows 1σ uncertainty of the stacked spectrum.

In the text
thumbnail Fig. 10.

Relation of EW(CIII]) with (from left to right) stellar mass, Ks-band luminosity, and FUV luminosity, for the stacks A, B, C, and D for each parameter, and for the individual galaxies in the C3 sample (small circles; in green, we show galaxies at z ≲ 3, and in red, galaxies at z ≳ 3). In the x-axis, the parameters correspond to mean values for stacks, while for the C3 sample, the parameters correspond to the SED fitting values.

In the text
thumbnail Fig. 11.

Diagnostic diagrams with EW of UV emission lines for the different stack sets, color-coded by the parameter used for stacking. The black dashed lines separate AGNs and star-forming galaxies as proposed in Nakajima et al. (2018a). The dotted and dotted-dashed lines are from Hirschmann et al. (2019), meant to separate between SFGs and composite as well as composite and AGN, respectively. Top row: diagnostic of EW(CIII])–CIII]/HeII ratio. The gray shaded region is where the models overlap. Middle row: diagnosis of EW(CIV)–CIV/HeII ratio. Bottom row: diagnosis of EW(OIII])–OIII]/HeII ratio. The green boxes are the composites from Le Fèvre et al. (2019) and the higher the size, the higher the EW (CIII]). The red X mark is the composite from Amorín et al. (2017). The red rectangles are the stacks from Nakajima et al. (2018b) with the brightest MUV (smallest rectangle), smaller EW(Lyα), faintest MUV, and larger EW(Lyα) (largest rectangle).

In the text
thumbnail Fig. 12.

Diagnostic diagrams with fluxes of UV emission lines for the different stack sets, color-coded by the parameter used for stacking. Symbols are the same as in Fig. 11. The red and cyan points are the photoionization models for AGN (Feltre et al. 2016) and SF (Gutkin et al. 2016) for metallicity Z = 0.0001, 0.0002, 0.0005, 0.001, 0.002, and for ionization parameter log(U) = −4, −3, −2, −1. The dotted-dashed lines in the bottom panel separates composite and AGNs according to Hirschmann et al. (2019).

In the text
thumbnail Fig. 13.

EW(Lyα)–EW(CIII]) relation. The results from stacks B are represented by the triangle symbols, with different color depending on the parameter used for stacking: blue for stellar mass, yellow for FUV luminosity, magenta for Ks luminosity, red for EW(CIII]), and green for EW(Lyα). Individual galaxies with measured Lyα from the C3 sample are the small blue circles with their typical errors on the upper left. Previous results from literature at similar redshift are also displayed from stacking (Shapley et al. 2003; Amorín et al. 2017; Nakajima et al. 2018b; Cullen et al. 2020; Feltre et al. 2020) and individual objects (Stark et al. 2014). Dashed black line is the best fit in Eq. (1) including our stacks and the sample from literature. Dashed red line is the best fit in Eq. (2), including just our stacks.

In the text
thumbnail Fig. 14.

EW(CIII])-stellar metallicity relation. The results from our A stacks are represented by star symbols. Colors for our stacks are the same as in Fig. 13. Black dashed line is the best fit in Eq. (3) and the gray shaded region is the 3-σ band uncertainty.

In the text
thumbnail Fig. 15.

EW(Lyα)-stellar metallicity relation. The symbols and colors for our stacks are the same as in Fig. 13. Black dashed line is the best fit to the stacked data (triangles) presented in Eq. (4), with the gray shaded regions at 3-σ, respectively. The black symbols are the stacks in Cullen et al. (2020) with Lyα in emission.

In the text
thumbnail Fig. 16.

Relation between the EW(CIII]) and C/O ratio as derived using HCM-UV. The red dashed line corresponds to the solar value. The cyan shaded region is the average C/O for Lyman-break galaxies in Shapley et al. (2003) and the red circle is the result for a composite spectrum of CIII]-emitters in Amorín et al. (2017). Individual CIII] emitters at high-z (red squares, Amorín et al. 2017) and low-z (black squares and circles, Berg et al. 2019a; Senchyna et al. 2021) are also included. Colors for our stacks are the same as in Fig. 13. The black dashed line is the best linear fit to our stacks, as displayed in Eq. (5).

In the text
thumbnail Fig. 17.

Relation between the different physical properties used for the stacks and C/O ratio using as derived using HCM-UV. The red dashed line corresponds to the solar value. The cyan shaded region is from Shapley et al. (2003). The results of the stacks are color coded by EW(CIII]). The black dashed line is the best linear fit for each parameter in Eqs. (6)–(8), respectively.

In the text
thumbnail Fig. 18.

Relation between SFR–EW(CIII]) (top panel) and SFR–C/O (bottom panel). Stars symbols are Stack A color-coded depending on the parameter used for stacking, according to legend, as in Fig. 13. Top panel: the black squares are results from CIII] emitters in Maseda et al. (2017). The small dots are the C3 sample with same colors as in Fig. 10. Bottom panel: the dashed black line is the best fit for the SFR–CO relation in Eq. (9).

In the text
thumbnail Fig. 19.

Relations between stellar mass and metallicity. Left panel: MZRs with the stacks A, B, D by stellar mass with star, triangle and circle symbols, respectively, and color-coded by EW(CIII]). The black symbols are stacks from Cullen et al. (2019, 2021), Calabrò et al. (2021). The orange circles is the local MZRs from Zahid et al. (2017). The blue dashed line is the best fit with the shaded region at 3-σ in Eq. (10), and the black dashed line is the MZRs from Cullen et al. (2019). Right panel: MZRg with the same stacks as left panel but re-scaled with the α-enhancement. Green dashed line from Steidel et al. (2014). Red dashed line from Troncoso et al. (2014). Black symbols are values from literature (Sanders et al. 2021; Cullen et al. 2021). The orange points are the local MZRg from Andrews & Martini (2013). The Zg scale in the y-axes corresponds to the gas-phase metallicity obtained from Eq. (11) (see text for more details).

In the text
thumbnail Fig. 20.

FMR with the stacks A, B, D by stellar mass, color-coded by EW(CIII]). Blue dashed line is the relation shown in Mannucci et al. (2010). The right scale is the [O/H] ∼ log(Zg/Z) assuming the α-enhanced in Eq. (11).

In the text
thumbnail Fig. 21.

Relation between C/O and metallicity. Left panel: C/O–Z relation. Black symbols are values of stars from the Galactic thin, thick disk, halo and unclassified from the literature (Amarsi et al. 2019). The dashed black line is the K20 model in Kobayashi et al. (2020) and the dotted line is the model MWG-11 in Romano et al. (2020). Right panel: C/O–Zg relation. The red markers are high-redshift galaxies from literature and black markers are local galaxies and HII regions (see text for references). The multi-zone chemical evolution models from Mattsson (2010) are also shown by black lines. The green dotted-dashed line is the chemical model from Mollá et al. (2015). Both panels: the red dashed line is the C/O solar value. The stack A by EW(CIII]) and broad band luminosities are shown by markers according to legend and color-coded by EW(CIII]).

In the text
thumbnail Fig. A.1.

MFUV–EW(CIII]) relation for the stacks by FUV luminosity with the non-detected CIII] sample (in cyan symbols). Faint symbols are the same as in Fig. 10.

In the text
thumbnail Fig. A.2.

EW(CIII])–C/O relation for the stacks by FUV luminosity with the non-detected CIII] sample (in cyan symbols). Faint symbols, shaded regions, and dashed lines are the same as in Fig. 16

In the text
thumbnail Fig. B.1.

Resulting Stack B for each physical parameter. Left to right: MKs, MFUV, and stellar mass. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line in the 1-σ error spectrum. The vertical lines mark known UV lines. Information about the number of galaxies, the mean redshift, and the mean parameter are included in each panel.

In the text
thumbnail Fig. B.2.

Resulting Stack C for each physical parameter. Left to right: MKs, MFUV, and stellar mass. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line in the 1-σ error spectrum. The vertical lines mark known UV lines. Information about the number of galaxies, the mean redshift and the mean parameter are included in each panel.

In the text
thumbnail Fig. B.3.

Resulting Stack D for each physical parameter. Left to right: MKs, MFUV, and stellar mass. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line in the 1-σ error spectrum. The vertical lines mark known UV lines. Information about the number of galaxies, the mean redshift and the mean parameter are included in each panel.

In the text
thumbnail Fig. B.4.

Stacks by EW(CIII]). Left to the right on top row: Stack A, Stack D. Middle row: Stack B. Bottom row: Stack C. In each panel, the green faint line is the stack spectrum with the ∼0.6 Å/pixel sampling, while the black one is with ∼1.2 Å/pixel. The blue line in the 1-σ error spectrum. The vertical lines mark known UV lines. Information about the number of galaxies, the mean redshift and the mean parameter are included in each panel.

In the text
thumbnail Fig. C.1.

Relation between EW(CIII]) and stellar mass (left panel) and FUV (right panel) luminosity, similar as shown in Fig. 10 for stack A (blue stars) and C3 sample (small blue circles). The cyan pentagons are the median values for the individual galaxies in the bin and their errorbars correspond to the standard deviation in the bin (boundaries marked by vertical dashed line).

In the text
thumbnail Fig. C.2.

Changes in linear fit in relations presented in Figs. 14, 15, 16, and 18, considering only stacks based on luminosities. Symbols are the same as in Figs. 14, 15, 16, and 18, respectively. The cyan dashed line shows the best fit including only the A stacks binned by luminosity. The stacks excluded in the fit are shown using transparent symbols.

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.