NOEMA reveals the true nature of luminous red JWST z > 10 galaxy candidates

The ﬁrst year of JWST has revealed a surprisingly large number of luminous galaxy candidates beyond z > 10. While some galaxies have already been spectroscopically conﬁrmed, there is mounting evidence that a subsample of the candidates with particularly red inferred UV colours are, in fact, lower redshift contaminants. These interlopers are often found to be ‘HST-dark’ or ‘optically faint’ galaxies at z ∼ 2 − 6, a population that is key to improving our understanding of dust-obscured star formation throughout cosmic time. This paper demonstrates the complementarity of ground-based mm-interferometry and JWST infrared imaging to unveil the true nature of red 1.5–2.0 µ m dropouts that have been selected as ultra-high-redshift galaxy candidates. We present NOEMA Polyﬁx follow-up observations of four JWST red 1.5–2.0 µ m dropouts selected by Yan et al. (ApJ, 942, L8) as ultra-high-redshift candidates in the PEARLS-IDF ﬁeld. The new NOEMA observations constrain the rest-frame far-infrared continuum emission and e ﬃ ciently discriminate between intermediate-and high-redshift solutions. We report > 10 σ NOEMA continuum detections of all our target galaxies at observed frequencies of ν = 236 and 252GHz, with FIR slopes indicating a redshift of z < 5. We modelled their optical-to-FIR spectral energy distribution (SED) with multiple SED codes, ﬁnding that they are not z > 10 galaxies but dust-obscured, massive star-forming galaxies at z ∼ 2 − 4 instead. The contribution to the cosmic star formation rate density (CSFRD) of such sources is not negligible at z (cid:39) 3 . 5 ( φ (cid:38) (1 . 9 − 4 . 4) × 10 − 3 cMpc − 3 ; or > 3 − 6% of the total CSFRD), in line with previous studies of optically faint and sub-millimeter galaxies. This work showcases a new way to select intermediate-to high-redshift dust-obscured galaxies in JWST ﬁelds with minimal wavelength coverage. This approach opens up a new window onto obscured star formation at intermediate redshifts, whilst removing contaminants with red colours from searches at ultra-high redshifts.


Introduction
The first year of operation of JWST has transformed our view of early galaxy evolution.In particular, the imaging sensitivity and near-infrared coverage up to ∼5 µm of the NIRCam instrument (Rieke et al. 2023) has been unveiling new galaxies undetected even in the deepest HST observations.These objects fall in two broad categories: 1) ultra-high-redshift candidates (z > 10) invisible in HST due to the Lyman-α break caused by the neutral intergalactic medium and 2) a collection of massive, intermediate-redshift (2 z 7) galaxies with red colours and high dust attenuation making them extremely faint at the wavelengths probed by HST ( 1.5 µm).
Ultra-high-redshift galaxies have been the focus of multiple studies in the first deep fields observed with JWST, with numerous candidates claimed at z > 10, as well as some up to z > 15 (Finkelstein et al. 2022(Finkelstein et al. , 2023a;;Naidu et al. 2022b;Harikane et al. 2023b;Donnan et al. 2022;Adams et al. 2022;Yan et al. 2023a,b).An excess at the bright end of the z 10 UV luminosity function could indicate enhanced star formation rate efficiencies in the first 500 Myr, but caution is still necessary as most of the objects still only have photometric redshifts.Whilst multiple objects have been successfully confirmed with spec-troscopy (e.g.Arrabal Haro et al. 2023a,b;Wang et al. 2023;Harikane et al. 2023a;Bunker et al. 2023a;Curtis-Lake et al. 2023), catastrophic outliers have already been identified, with z ∼ 5 galaxies with extreme dust and/or line emission properties masquerading as z > 10 candidates (Naidu et al. 2022a;Zavala et al. 2023;Arrabal Haro et al. 2023a;Wang et al. 2023).Some studies, recognising the likely contamination of highredshift samples by lower redshift dusty galaxies, have already proposed additional criteria (e.g.∆χ 2 between low-and highz solutions, integral of p(z > 7, 8, 9)) to complement the initial photometric redshift selection (e.g.Finkelstein et al. 2023a;Harikane et al. 2023b).These different definitions of robust high-redshift candidate lead to disagreements on the exact number of candidates in each field (Bouwens et al. 2023).Whilst such ad hoc prescription have successfully led to large fraction of objects confirmed with spectroscopy (e.g.Finkelstein et al. 2023b;Harikane et al. 2023a), the nature of the dusty, likely lower-redshift contaminants remains speculative in the absence of dedicated follow-up observations.The JWST-detected population of dusty 3 z 6 galaxies has been identified as an extension of 'HST-dark', 'ALMA-only', optically faint infrared, or dusty star-forming galaxies (DSFGs) thought to be responsible for most of the obscured star formation at high redshift (e.g.Blain et al. 1999;Smail et al. 2002;Chapman et al. 2005;Barger et al. 2012;Swinbank et al. 2014;Elbaz et al. 2011;Casey et al. 2014Casey et al. , 2018;;Wang et al. 2019;Alcalde Pampliega et al. 2019;Williams et al. 2019;Fudamoto et al. 2020;Dudzevičiūtė et al. 2020;Shu et al. 2022;Xiao et al. 2023).Their detection and study with JWST shows that these names only described a particular selection function for objects with a wide range of properties (Nelson et al. 2023;Rodighiero et al. 2022;Pérez-González et al. 2023;Barrufet et al. 2023;Smail et al. 2023;Barger & Cowie 2023).Broadly speaking, these studies find that optically faint galaxies detected now with JWST are intermediate-to high-redshift (2 z 8), dusty (A v > 2) with a large fraction of obscured star formation or passive galaxies, and could contribute a significant fraction >10% of the cosmic star formation rate density at z 5.However, this picture is evolving rapidly as their census is still incomplete and new ways of identifying such optically invisible galaxy are still in development.
In particular, large numbers of these objects might be hidden in ultra-high-redshift candidates lists, based on the detection of a break at 1−2 µm.The purpose of this Letter is to demonstrate the power of millimeter observations to determine the nature of ambiguous 1.5-2.0µm red dropouts detected in JWST imaging data.Specifically, we present follow-up observations of four 1.5-2.0µm red dropouts selected by Yan et al. (2023a) in the PEARLS survey with a limited number of JWST/NIRCam imaging filters (F150W, F277W, F356W, and F444W).We show that most SED fitting codes, using standard parameters, either prefer z > 10 or have degenerate or bimodel posterior redshift solutions.Here, we present NOEMA 1 mm detections of their FIR continuum (100% detection rate) to demonstrate that such observations break the degeneracy between intermediate-and high-redshift photometric redshift solutions.Our SED modelling of the red 1.5-2.0µm dropouts shows that they are massive, dusty galaxies at 2 < z < 4, similar to so-called 'HST-dark' galaxies (see also Zavala et al. 2023).

NOEMA observations
We studied four objects detected in JWST F200W, F356W, and F444W imaging of the PEARLS survey (Windhorst et al. 2023) and selected as ultra-high-redshift candidates based on their colours (Yan et al. 2023a).The four targets of this work were chosen to be the brightest, but also the reddest in the Yan et al. (2023a) sample (F150D_JWIDF_E01,F150D_JWIDF_H17, F200D_JWIDF_M03, and F200D_JWIDF_H08), guaranteeing a S /N ∼ 10 detection at ∼250 GHz if their redshift is z ∼ 3−5, or S /N >∼ 100 if they are at z > 10 (see further Sect.3 and Appendix A).For completeness, their HST and JWST photometry from Yan et al. (2023a) is reproduced in Table 1, alongside our NOEMA FIR continuum measurements.
The targets were observed with NOEMA in Band 3 between December 17-18 2022.The targets were set at the phase center and observed for 0.5h (on-source) in track-sharing mode.Two targets are close (potentially interacting/merging; F150DE01 and F200DH08) and are observed simultaneously in the same pointing.The spectral setups were tuned at a nominal central frequency of 251.76 GHz in the upper sideband and 236.25 GHz in the lower sideband.
The data was calibrated and reduced at IRAM (remotely) using the latest CLIC package in the GILDAS framework.We flag a small fraction of the tracks affected by poor weather conditions and baselines with bad phase solutions.We image the continuum in the lower and upper sidebands separately using natural weighting and Hogbom cleaning down to 2σ in MAPPING with r = 5 circular support regions centered on the JWST target positions.The final synthesized beam size is 0.96 × 0.83 (0.90 × 0.77 ) and the final continuum rms achieved is 62−64(73−76) µJy beam −1 at 236.25(251.76)GHz, respectively.
We show the continuum imaging of the targets in Fig. 1, where all the sources are clearly detected at S /N 10-15.No additional continuum sources were detected in the NOEMA field of view (θ ∼ 30 ).The sources are only marginally resolved and the difference between the aperture-integrated fluxes and that of the central pixel is at most 10−20%, depending on the aperture size used.Thus, we simply measured the continuum fluxes in the central brightest pixel, which is equivalent to the total flux within one (synthesised) beam.This also avoids any complications when comparing the fluxes between the two sidebands.
The two close galaxies F150DE01 and F200DH08 are not resolved, with a total continuum flux of f 236 GHz = 1.04 ± 0.06 mJy, f 252 GHz = 1.19 ± 0.07 mJy.Based on the assumption that both objects are of a similar nature, we assign half of the flux to each galaxy and use a fiducial relative error of 34% reflecting the uncertainty on the fraction of flux originating from each object.This has little impact on our results as the redshift and nature of the sources is already sufficiently constrained by 10 mJy continuum detection.The continuum flux densities and errors for all sources are presented in Table 1.
We also produced datacubes for each sideband with 20 MHzwide channels.The rms noise in the cube 1.2(1.4)mJy beam −1 per channel.For each source, we extract spectra in r = 2 apertures as well as in the central pixel only and we find no emission line in the frequency range covered by our observations (232.5−240GHz, 248−255.5 GHz).(Boquien et al. 2019), BAGPIPES (Carnall et al. 2018), and PROSPECTOR (Johnson et al. 2021).We briefly review the parameters used to run each code below.

Multi
EAZY.The most important choice when using EAZY is the choice of templates.We ran EAZY twice: once with the tweak_FSPS_templates, augmented with the three bluer templates (with and without Lyman-α emission) from Larson et al. (2023), and a second time with the blue_sfhz1 templates with redshift-dependent star formation histories.We fit the data between redshift 0 and 20, using no prior on the luminosity, and applied template error using the values derived from the COSMOS2020 dataset.
CIGALE.We ran a grid of CIGALE models using a delayed star formation rate history (0.1 < τ/[Gyr] < 5) and a recent burst, using the SSP models from Bruzual & Charlot (2003) with a Chabrier (2003) IMF and a fixed metallicity of Z = 0.02.We varied the nebular line emission contribution with −4 < log U < −1 (in steps of 0.5) with a fixed gas metallicity of Z gas = 0.02 and fixed electron density of n e = 100 cm −3 .The dust attenuation was modelled with the dustatt_modified_starburst module that uses the Calzetti et al. (2000) attenuation for the continuum, with 0 < A V < 8.2 and a MW-like attenuation curve for the emission lines.The models were fitted within 0.0 < z < 20 in steps of ∆z = 0.05.
We also fit the near-infrared to FIR photometry with CIGALE templates, including an AGN contribution using SKIRTOR templates (Stalevski et al. 2016).We find that the addition of AGN template does not improve the χ 2 of the best-fit solutions and that the contribution of AGN to the observed photometry is negligible in the best-fit composite SEDs.We thus proceed by fitting the available photometry with pure stellar light SEDs.
BAGPIPES.We used the latest version of BAGPIPES (Carnall et al. 2018) with a delayed star formation history (0.1 < τ/[Gyr] < 10), the dust attenuation from Calzetti et al. (2000), with 0 < A v < 8, nebular line emission with −4 < log U < −2, metallicities spanning 0−2.5 Z .BAGPIPES uses Bruzual & Charlot (2003) SSP models with a Kroupa (2001) IMF, and we use a fixed metallicity Z = 0.02, matching that used in CIGALE.As BAGPIPES does not include a special treatment of non-detections, we set the fluxes of upper limits to 0 and use the 1σ upper limit as the error when fitting the spectra.
Prospector.We ran Prospector using a uniform redshift prior 0.1 < z < 20 and a delayed-tau SFH history.The physical parameters priors follow the choices of Tacchella et al. (2023), for instance the gas and stellar metallicity, dust properties, and nebular emission parameters (see Table 1 and Sect.3.4 of Tacchella et al. 2023), except for the dust extinction A V , where we used a uniform prior 0 < A V < 8 to match that used for the BAGPIPES and CIGALE runs above.
We fit the photometry of each object (see Table 1) with and without our NOEMA 1 mm constraints.We show in Fig. 2 the best-fit SEDs and the posterior redshift distributions.The first result is that a majority of the codes, in the absence of FIR constraints, either prefer a high-redshift (z > 6) solution or allow one (with a z > 10 solution peaking at least at half the likelihood of the intermediate redshift solution).The notable exception is EAZY using the new blue_shfz_13 templates, which consistently prefers the lower redshift solution, although L3, page 3 of 8 Fig. 2. SED and posterior redshift distributions of the four red 1.5-2.0µm dropouts studied in this work.The best high-redshift SED using only the JWST+HST photometry are shown with dotted lines and the best-fit low-redshift solutions found including the NOEMA data are shown with full lines (left).In the absence of FIR constraints, the four JWST bands and the one HST (black squares and limits) band cannot constrain the redshift distribution and allow or prefer a z > 10 for most sources and codes (middle).The addition of the FIR constraints (red circles) break the degeneracies and clearly reveal the z ∼ 2−4 nature of these sources (right).Upper limits are shown at the two σ level.
L3, page 4 of 8 2. Key physical properties of the four 1.5-2.0µm red dropouts.Notes.The properties are derived using the full optical to FIR SED constraints.For each parameter we give the median value and error from the 16th, 50th, and 84th percentiles of the BAGPIPES posterior and CIGALE's Bayesian estimates.
the posterior redshift distribution for F150W dropouts still contains solutions at z > 10.The different performance of the blue_shfz_13 templates can be explained by the inferred restframe UV slopes of our targets if they would be at z > 10.To derive the inferred rest-frame UV β slopes, we fit the rest-frame UV (1265−2580 Å) of all the z > 10 models from BAGPIPES, with a power law of f λ ∝ λ β .We find UV slopes for the z > 10 solutions β = 0.13 +0.11 −0.13 , −0.12 +0.17 −0.16 , 1.28 +0.21 −0.36 , 1.44 +0.12 −0.12 , for F150D_JWIDF_E01, F150D_JWIDF_H17, F200D_ JWIDF_H08, and F200D_JWIDF_M03, respectively.These β values are much higher than that the typical −3 β −1 measured in confirmed z > 10 galaxies (e.g.Bunker et al. 2023b;Curtis-Lake et al. 2023).The absence of extremely red templates at high-redshift in the blue_shfz_13 set therefore explains the lower photometric redshift solutions for this particular EAZY run.
Prospector also prefers lower-redshift solutions, although the probability for a z > 10 is still significant (p(z > 10) ∼ 0.1−0.4,see Fig. 2).However, we note that the Prospector results heavily depend on the choice of priors for the various physical parameters.Indeed, the maximum-likelihood solutions (plotted in dotted purple in Fig. 2) are z ML = 10.9, 3.66, 13.8, 15.6 for F150D_JWIDF_E01, F150D_JWIDF_H17, F200D_JWIDF_H08, and F200D_JWIDF_M03, respectively.We also find the posterior redshift distribution to be strongly dependent on the dust parameters priors, for instance, a lognormal prior for the dust extinction will result in z > 10 solution being strongly preferred.The Prospector results are likely priordominated due to the low number of continuum datapoints.
The second important result is that the addition of the NOEMA observations completely transforms the posterior redshift distribution, as also shown in Zavala et al. (2023).Indeed, once the FIR continuum constraints are included in the fits, an intermediate redshift solution (z ∼ 2−4) is strongly preferred in all codes2 .
Using BAGPIPES and CIGALE to fit their optical-to-FIR SED, we find that the four sources are galaxies at intermediate redshift (z ∼ 2−4), with star formation rates of SFR ∼ Fig. 3. Inferred star formation rate against stellar mass for the galaxies studied in this work.Results from BAGPIPES are shown in black circles, whereas those from CIGALE are shown in grey squares.The SFR upper limit (F200D_JWIDF_M03 fitted with CIGALE) is shown at the 2σ level.The coloured lines show the galaxy main sequence at z = 2, 3, 4 from Schreiber et al. (2015).20−150 M yr −1 , stellar masses of 9.7 < log M * /[M ] < 11.4, and high obscuration (2 A V 7).CIGALE tends to prefer solutions with lower redshifts, along with higher masses and dust attenuations, but the BAGPIPES and CIGALE results are consistent within the 1σ errors (see further Table 2 for the full results of the two codes).The most important conclusion is that these objects are obscured galaxies at z ∼ 2−4 that lie on the galaxy main sequence (see Fig. 3), similarly to the HST-faint/JWSTdetected objects reported by Barrufet et al. (2023).
We used the best-fit BAGPIPES spectrum to derive the obscured and unobscured star formation rate, using the Kennicutt & Evans (2012) conversion factors from the restframe FIR (3−1100 µm) and UV (1550 Å) luminosities.We find that that all galaxies are highly obscured, with a negligible contribution of the unobscured SFR to the total SFR, with SFR UV /(SFR UV + SFR IR ) ∼ 10 −7 −10 −3 .We also estimate the contribution of these objects to the total cosmic star formation rate density using the median SFR of the BAGPIPES/CIGALE posterior and the area of the PEARLS-IDF field our targets were selected from (14.2 arcmin 2 ).Using a wide redshift bin ∆z = 1 and ignoring the effect of the unknown selection function, we find lower limits on the cosmic SFR density (CSFRD) of φ = (3.2−13.0)× 10 −4 M yr −1 cMpc −3 at z = 2 and φ = (1.9−4.4) × 10 −3 M yr −1 cMpc −3 at z = 3.5, where the range of values encompasses the scatter between the CIGALE and BAGPIPES results (see Table 2).The CSFRD we derive for the four sources is in good agreement with the literature constraints on the contribution of obscured galaxies to the CSFRD at these redshifts (see Fig. 4, right panel;and e.g. Blain et al. 1999;Chapman et al. 2005;Barger et al. 2012;Swinbank et al. 2014;Wang et al. 2019;Dudzevičiūtė et al. 2020;Shu et al. 2022;Barrufet et al. 2023), representing a >0.2−1.0(>3−6%)contribution to the total CSFRD at z = 2(3.5).

Conclusions
We present a NOEMA follow-up study of red 1.5-2.0µm dropouts in JWST PEARLS-IDF, showing that sub-mm L3, page 5 of 8 Fig. 4. Cosmic star formation rate density contribution from dusty starforming galaxies.The constraints from this work sample are shown with red squares (BAGPIPES) and circles (CIGALE).Selected constraints from the litterature are shown in black (Swinbank et al. 2014;Wang et al. 2019;Barrufet et al. 2023, but see further text).We also show the total SFRD with a dashed black line (Madau & Dickinson 2014).The extremely red JWST dropouts represent 5% of the CSFRD at z 3−4 follow-up is highly efficient at unveiling their true nature.Their FIR continuum is detected at a signal-to-noise ration of S /N ∼ 15 in our the 236, 252 GHz observations, with typical values of 0.5−1.5 mJy.Our SED modelling shows that such a faint continuum in the mm-regime unambiguously argues against an ultra-high-redshift (z > 10) nature (Fig. A.1 and Appendix A).Instead, the 1.5-2.0µm breaks and red JWST colours are shown to select massive, star-forming, and dust-obscured galaxies at z ∼ 2−4, extending the main sequence of SMGs and optically faint galaxies.These objects contribute a non-negligible fraction of the CSFRD at z 3.5 (φ (1.9−4.4)×10−3 M yr −1 cMpc −3 ), as found in previous studies specifically targeting optically faint 'HST-dark' or DSFGs galaxies.This demonstrates the efficiency of combining mm continuum observations with infrared imaging in a limited number of JWST filters to probe obscured starformation beyond cosmic noon.
This work show specific cases where a simple ∼1−2 µm break selection is prone to select optically faint galaxies instead of ultra-high-redshift galaxies when no further criteria on the UV slopes of star-formation rate histories have been applied.Fortunately, we have demonstrated that red objects can be robustly identified with observations of their rest-frame FIR (see also Zavala et al. 2023).The expected peak FIR continuum of ∼1 mJy (and the characteristic continuum slope; see Fig . A.1) can be detected at 10σ in ∼3(20) min of on-source time with ALMA (NOEMA), representing a highly efficient alternative to spectroscopic follow-up with JWST.We note that this is a different situation to previous works focusing on JWST high-redshift objects with blue colours.There, expensive redshift looking for rest-frame FIR emission lines must be conducted as the continuum is likely too faint to be detected and does not significantly constrain the redshift of the galaxy (e.g.Kaasinen et al. 2023;Bakx et al. 2023;Popping 2023;Fujimoto et al. 2023).
We have also demonstrated that the use of bluer templates at higher redshift, such as those made recently available with EAZY (blue_sfhz_13), improve the photo-z selection and characterisation of red dropouts in the absence of observations in the mm regime.However, we note that these templates find a nonnegligible high-redshift probability for the two F150W dropouts (e.g.p(z > 10) = 0.18 and 0.35, see Fig. 2).Similarly, Prospector is found to prefer lower-redshift galaxies even in the absence of FIR constraints.Yet it also gives broad posterior redshift distributions, which allow for a z > 10 solution, even when including our NOEMA constraints.Whilst some studies have already rejected such objects with a substantial low-redshift probability from z 10 samples with a variety of additional criteria (e.g.Finkelstein et al. 2023a,b;Harikane et al. 2023b), further tests of these templates and redshift priors with statistical samples of spectroscopically-confirmed red JWST-only galaxies are still necessary to improve the robustness of photometric redshifts.Another unexplored avenue is the use of morphology priors, as in the case of the red objects studied in this work or those in Nelson et al. (2023), for instance, which have elongated or disky morphologies; these are evidently in stark contrast with the almost point-like z > 10 sources confirmed so far.Finally, for surveys with a large enough field of view, a CSFRD prior could be imposed to guide the ensemble prediction for all galaxies in the field (see more in Appendix A).
In summary, distinguishing obscured star-forming galaxies from ultra-high-redshift galaxies with sparse JWST wavelength coverage is difficult when using standard SED codes (templates), although dedicated templates and/or appropriate priors could improve some of the photometric redshifts.However, follow-up (sub-)mm observations of the FIR continuum provide a highly efficient way to determine the nature of high-redshift dusty starforming galaxies.
-wavelength SED modelling and the nature of red 1.5-2.0µm dropouts We modelled the SED of the four red dropouts using five different codes: EAZY (Brammer et al. 2012), CIGALE L3, page 2 of 8

Fig. 1 .
Fig. 1.NOEMA continuum observations of the four red 1.5-2.0µm dropouts from the PEARLS-IDF data in Yan et al. (2023a, indicated with cyan symbols).We show the 236 GHz (top row) and 252 GHz continuum (bottom row) for each object with the galaxies position from JWST indicated cyan symbols.The full contours indicate the (3, 6, 9, 15)σ level, where σ is the root mean square (rms) noise in the images.The dashed contours show the −3σ level.The red dashed circle shows a r = 2 aperture.

Table 1 .
Yan et al. (2023a)y of the PEARLS-IDF targets studied in this work.The measured HST and JWST photometry is reproduced fromYan et al. (2023a)for completeness.Upper limits are given at the 2σ level.